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1.  Formulation  and  Applications 

An  economic  problem  of  great  practical  importance  is  to  choose  the 
location  of  facilities,  such  as  industrial  plants  or  warehouses,  in  order  to 
minimize  the  cost  (or  maximize  the  profit)  of  satisfying  the  demand  for  some 
commodity.  In  general  there  are  fixed  costs  for  locating  the  facilities  and 
transportation  costs  for  distributing  the  commodities  between  the  facilities 
and  the  clients.  This  problem  has  been  extensively  studied  in  the  literature 


and  i s  often  referred  to  as  thej?1ant. ^warehouse  or  facility 


location 


problem.  When  each  potential  facility  has  a  capacity,  which  is  the  maximum 

"r 

demand  that  it  can  supply,  the  problem  is  known  as  the  capacitated  facility 

V 

location  problem.  When  the  capacity  hypothesis  is  not  needed,  we  have  the 
simple  or  uncaoacitated  facility  location  problem,  which  w6  abbreviate  by  UL. 

P*  -*■ 

The  mathenatical  formulation  of  these  problems  as  integer  programs  has 

proven  very  fruitful  in  the  derivation  of  solution  methods.  To  formalize  UL, 

r\ 

we  consider  a  set  of  clients  I  »  {l,...,m>  with  a  given  demand' fpr  a  single 
commodity,  and  a  set  of  sites  J  *  {I,...,n}  where  facilities  can  be  - 
located.  In  the  literature  it  has  been  traditional  to  use  the  phrasing  x  , 
"facility  j  is  open"  to  mean  that  a  facility  is  actually  located  at  site 


\ 


j.  Let  f,  be  the  given  fixed  cost  of  opening  facility  j  and  assume  there 

J 

is  a  known  profit  c ■ .  that  is  made  by  satisfying  the  demand  of  client  i 

*  J 

from  facility  j.  Typically,  c-j  is  a  function  of  the  production  costs  at 
facility  j,  the  transportation  costs  from  facility  j  to  client  i,  the 
demand  of  client  i  and  the  selling  price  to  client  i.  For  example, 
cij  *  d i ( P, -t 1 j )  where  d^  is  the  demand,  p.  the  price  per  unit,  qj 
the  production  cost  per  unit  and  t^.  the  transportation  cost  per  unit.  UL 


is  to  open  a  subset  of  facilities  in  order  to  maximize  total  profit,  given 

that  all  demand  has  to  be  satisfied. 

For  any  given  set  S  of  open  facilities,  it  is  optimal  to  serve  client 

i  from  a  facility  j  for  which  c^  is  maximum  over  j  «  S.  So,  given  S, 

the  profit  is  z(S)  *  >  max  -  V  f4.  The  problem  is  to  find  a  set  S 

id  j«S  jcS  J 

that  yields  the  maximun  profit  Z,  i.e.  Z  ■  max  z(S).  This  can  be  viewed  as 

ScJ 

a  combinatorial  formulation  of  the  problem.  Note  that  there  is  no  loss  of 
generality  in  assuming  f j  >_  0  for  all  j  e  J  since  if  ffc  <  0  every 
optimal  solution  contains  facility  k. 

An  integer  linear  programming  formulation  is  obtained  by  introducing  the 
following  variables.  Let  x^  *  1  if  facility  j  is  open  and  x^  *  0 
otherwise;  y^  *  1  if  the  demand  of  client  i  is  satisfied  from  facility 
j  and  y.j j  »  0  otherwise.  The  integer  program  is 


(1.1) 

Z  *  max  l  l  c^y^ 
iel  jeJ 

• 

1 

l  Vi 

jeJ  J  J 

(1.2) 

I  yi , 

je  J 

[  * 

1 

all  iel 

(1.3) 

yi; 

li 

Xj 

all  i  e  I,  j  e  J 

(1.4) 

xr  yu 

1 6 

(0,1)  all 

i  e  I,  j  e  J. 

The  constraints  (1.2)  guarantee  that  the  demand  of  every  client  is  satisfied 
whereas  (1.3)  guarantees  that  the  clients  are  supplied  only  from  open 


The  selling  price  of  the  commodity  to  client  1  enters  the  solution  of 
(1.1)-(1.4)  as  a  constant.  Only  the  costs  (production,  transportation  and 
fixed  operating  costs)  are  relevant  for  the  decision.  This  Is  why  In  the 
literature  UL  Is  often  presented  as 


d.l’) 


min  l 


-  I  dnyii +  I  Vi 

1«I  jc  J  1J  j«J  3  3 


subject  to  (1.2) -(1.4)  where  d^.  Is  the  cost  (production  plus  transpor¬ 
tation)  of  serving  client  1  from  facility  j.  This  formulation  Is  math¬ 
ematically  equivalent  to  (1.1)-(1.4)  since  (1.1* )  becomes  -(1.1)  by  setting 
C1j  *  “^Ij*  *n  othep  words»  costs  can  simply  be  regarded  as  negative  profits. 

In  the  integer  program  (1.1)-(1.4),  the  x^'s  are  the  Important  vari¬ 
ables.  Once  a  set  of  Xj's  that  satisfy  (1.4)  are  specified,  it  is  simple 

to  determine  a  set  of  y^'s  that  solves  the  Integer  program  for  the  fixed 

Xj's.  Even  If  we  drop  the  Integrality  requirement  on  the  y^'s,  an  optimal 

set  of  y^j's  Is  given  by  y^#  ■  1  and  y^  »  0,  j  *  j*  where  c^*  * 

max{Cjj:  x^  ■  1,  j  «  J).  Thus  (1.4)  can  be  replaced  by 


(1.4* ) 


Xj  e  (0,1),  y^  >  0  all  i  e  I,  j  e  J. 


The  formulation  (1.1)-(1.3),  (1.4')  Is  a  mixed  integer  linear  program. 


An  integer  linear  program  equivalent  to  (1.1M1.4)  is  obtained  by 
replacing  the  constraints  (1.3)  by  the  more  compact  set  of  constraints 
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m 
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( 1 -3  * ) 


I  yH  <  mx,  all  j  e  J. 
iel  J  J 


To  see  that  (1.3)  can  be  replaced  by  (1.3'),  note  that  when  x.  *  0  both 

J 

(1.3)  and  (1.3* )  imply  y.^  *  0  for  all  i  e  I,  and  when  x^.  *  1  both  (1.3) 
and  (1.3*)  are  satisfied  for  all  choices  of  y^  that  satisfy  the  constraints 
(1.2).  However,  the  two  formulations  are  equivalent  only  for  0-1  values  of 
the  variables  x^..  For  each  j,  (1.3* )  is  obtained  by  sunming  (1.3)  over  all 
i  e  I;  hence  any  solution  to  (1.2),  (1.3)  is  also  a  solution  to  (1.2)  and 
(1.3*).  But  the  converse  is  false  when  0  <  Xj  <  1;  i.e.  the  feasible  region 
defined  by  (1.2),  (1.3* )  and 


(1.5) 


®  <  «j  <  i.  *ij  >  °  ( « i « J 


strictly  contains  the  region  defined  by  (1.2),  (1.3)  and  (1.5). 

In  UL,  the  nunber  of  facilities  that  are  open  in  an  optimal  solution 
is  not  specified;  it  is  determined  by  a  solution  of  (1.1)-(1.4).  From  a 
practitioner's  standpoint,  it  might  be  useful  to  consider  a  formulation  where 
the  number  p  of  open  facilities  is  a  parameter  of  the  problem.  This  is 
realized  by  adding  one  of  the  following  constraints  to  the  program  (l.l)-(l .4) 


where  p  Is  some  given  Integer,  1  £  p  £  n.  The  formulation  (1.1)-(1.4), 

(1.6)  has  been  called  the  p-facllity  location  problem.  When  n  »  m  and 
f j  *  0  for  all  j  (  J,  (1.1)-(1.4),  (1.6)  Is  known  as  the  p-medlan  problem. 

It  Is  often  solved  in  the  literature  using  the  minimization  objective  func¬ 
tion  (1.1')  where  (d^)  is  a  distance  matrix,  i.e.  it  is  positive,  symmetric 
and  djj  +  d^  >.  d^k,  1  ^  1  ,j,k  jC  n.  This  model  Is  developed  in  Chapter  3  of 
the  book. 

To  motivate  further  the  study  of  UL,  we  give  two  examples.  First  we 
Interpret  (1.1)-(1.4)  as  a  bank  account  location  problem,  see  Cornuejols, 
Fisher,  Nemhauser  (1977b).  The  second  example  occurs  In  clustering  analysis, 
see  Mulvey  and  Crowder  (1979).  Other  examples  arise  in  lock-box  location 
[Kraus,  Janssen  and  McAdams  (1970)3  location  of  offshore  drilling  platforms 
[Balas  and  Padberg  (1976)],  economic  lot  sizing  [Krarup  and  Slide  (1977)], 
machine  scheduling  and  Information  retrieval  [Hansen  and  Kaufman  (1972)]  and 
portfolio  management  [Beck  and  Mulvey  (1982)]. 

A  Sank  Account  Location  Problem 

The  number  of  days  required  to  clear  a  check  drawn  on  a  bank  In  city  j 
depends  on  the  city  1  in  which  the  check  Is  cashed.  Thus,  to  maximize  its 
available  funds  a  company  that  pays  bills  to  clients  in  various  locations  may 
find  it  advantageous  to  maintain  accounts  in  several  strategically  located 
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banks.  It  would  then  pay  bills  to  clients  in  city  i  from  a  bank  in  city  j 
that  had  the  largest  clearing  time.  The  economic  significance  to  large 
corporations  of  such  a  strategy  is  discussed  in  an  article  in  Businessweek 
(1974). 

To  formalize  the  problem  of  selecting  an  optimal  set  of  account 
locations,  let  I  *  {l,...,m}  be  the  set  of  client  locations,  J  3  (l,...,n) 
the  set  of  potential  account  locations,  f  •  the  fixed  cost  of  maintaining  an 

J 

account  in  city  j,  di  the  dollar  volune  of  checks  paid  in  city  i,  and 
the  number  of  days  (translated  into  monetary  value)  to  clear  a  check  issued 
in  city  j  and  cashed  in  city  i.  All  this  information  is  assuned  to  be 
known  and  c^  *  di*ij  represents  the  value  of  paying  clients  in  city  i 
from  an  account  in  city  j.  Let  x.  3  1  if  an  account  is  maintained  in  city 

J 

j,  and  Xj  3  0  otherwise;  y^  3  1  if  clients  in  city  i  are  paid  from 
account  j,  and  y^.  3  0  otherwise.  Then  the  account  location  problem  can 
be  stated  as  ( 1 .!>-( 1 .4) . 

Besides  desiring  to  delay  paynents  for  as  long  as  possible,  corporations 
also  want  to  collect  funds  due  to  them  as  quickly  as  possible.  This  can  be 
done  by  situating  check  collection  centers  or  "lock-boxes"  at  optimal 
locations.  Since  (1.1)  and  (1.1')  are  mathematically  equivalent,  (1.1)-(1.4) 
is  also  a  model  for  the  lock-box  problem. 

A  Clustering  Problem 

Cluster  analysis  consists  of  partitioning  objects  into  classes,  known  as 
clusters,  in  such  a  way  that  the  elements  within  a  cluster  have  a  high  degree 
of  natural  association  among  themselves  while  the  clusters  are  relatively 
distinct  from  one  another.  Cluster  analysis  is  used  in  biology,  psychology. 


■a 


medicine,  artificial  Intelligence,  pattern  recognition,  marketing  research, 
automatic  classification  and  statistics. 

Let  I  ■  {l,...,m>  be  the  set  of  objects  to  be  clustered.  The 
clustering  is  done  around  objects  that  “represent"  the  clusters.  Let 
J  «  {l,...,n>  be  the  set  of  eligible  “cluster  representatives".  In  many 
applications  J  *  I.  However  In  some  applications  |j|  <  |l|,  l.e.  only 
objects  with  certain  characteristics  can  qualify  as  representatives  (e.g. 
survey  papers  and  books  In  the  automatic  classification  of  technical  material 
in  some  field).  In  other  cases  |l|  <  |j|,  (e.g.  in  the  automatic  classi¬ 
fication  of  technical  material,  an  alternative  to  the  above  policy  Is  to 
represent  a  cluster  by  a  list  of  key  words.  This  list  does  not  need  to 
match  exactly  that  of  any  single  paper  In  the  cluster). 

The  clustering  problem  Is  defined  relative  to  a  matrix  of  parameters 
c^  that  represent  the  similarity  between  objects  1  and  j.  For  example 
Cjj  could  be  the  number  of  common  key  words  associated  with  technical  papers 
1  and  j.  Anderberg  (1973)  gives  several  ways  of  calculating  the  similar¬ 
ity  matrix  C.  In  most  applications  there  are  no  natural  fixed  charges. 
Rather,  the  constraint  (1.6)  Is  added  to  the  formulation  (1.1)-(1.4). 
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2.  Brief  Historical  Overview 

There  is  a  vast  literature  on  the  uncapacitated  location  problem. 

Krarup  and  Pruzan  (1983)  give  an  up-to-date  survey.  We  will  not  repeat  their 
effort  here.  However,  to  set  the  stage  for  the  ensuing  sections,  we  will 
mention  the  main  solution  approaches  and  cite  some  basic  references. 

The  first  approaches  to  solve  UL  were  heuristic.  One  of  the  earliest 
heuristics  Is  due  to  Kuehn  and  Hamburger  (1963),  who  actually  present  it  for 
a  wider  class  of  location  problems.  It  consists  of  two  routines.  The  first 
routine  opens  facilities  sequentially  in  an  order  that  maximizes  the  increase 
of  the  objective  function  at  each  step.  It  stops  when  adding  a  new  facility 
could  only  decrease  the  objective.  Kuehn  and  Hamburger  called  it  the  "add 
routine";  In  the  modern  literature  such  a  procedure  is  called  a  greedy 
heuristic  because  of  Its  appetite  for  maximum  Improvement  at  each  step. 

Their  second  routine  is  the  "bump  and  shift  routine".  It  eliminates  (bumps) 
any  facility  that  has  become  uneconomical  because  of  the  presence  of  other 
facilities  chosen  subsequently  by  the  greedy  heuristic.  Then,  starting  from 
this  feasible  solution.  It  considers  interchanging  (shifting)  an  open  and 
closed  facility.  Such  a  pairwise  Interchange  Is  performed  if  It  Improves 
the  current  feasible  solution  and  the  procedure  stops  when  a  solution  has 
been  found  that  cannot  be  improved  by  such  interchanges.  In  the  remainder  of 
this  chapter,  the  shifting  procedure  will  be  referred  to  as  an  interchange 
heuri stic. 

The  greedy  and  interchange  heuristics  are  the  basis  of  numerous 
approximation  algorithms.  They  can,  of  course,  be  helpful  in  exact 
algorithms  that  require  feasible  Initial  solutions  or  use  feasible  solutions 
in  other  ways,  Spielberg  (1969b)  and  Hansen  and  Kaufman  (1972).  However, 


when  a  heuristic  Is  used  to  obtain  a  final  solution,  it  is  very  Important  to 
have  an  upper  bound  as  well  so  that  the  user  can  be  confident  that  the 
heuristic  solution  value  is  not  too  far  from  the  optimal  value.  Cornuejols, 
Fisher  and  Nemhauser  (1977b)  gave  such  bounds  for  the  greedy  and  interchange 
heuristics,  both  a  priori  worst-case  bounds  and  a  posteriori  bounds  on  a 
particular  solution  constructed  by  the  heuristic.  These  bounds  were 
generalized  to  the  maximization  of  submodular  set  functions  (see  Section  9) 
by  Nemhauser,  Wolsey  and  Fisher  (1978). 

General  solution  techniques  for  finding  optimal  solutions  to  integer 
programs  have  been  customized  for  UL.  The  mixed  integer  linear  programming 
formulation  can  be  solved  by  Benders  decomposition.  Benders  (1962).  This 
approach  was  proposed  by  Bal inski  and  Wolfe  (1963)  and  appears  to  have  been 
the  first  attempt  to  solve  UL  to  optimality.  The  computational  experiments 
were  discouraging,  see  Bal Inski  (1965),  and  this  method  was  abandoned  until 
recently.  Magnanti  and  Wong  (1977)  develop  techniques  to  accelerate  the 
convergence  of  Benders  decomposlton.  They  generate  strong  cuts  from  the  set 
of  feasible  Benders  cuts  and  by  so  doing  they  are  able  to  reduce  the  number 
of  integer  programs  to  be  solved.  Nemhauser  and  Wolsey  (1981)  consider  the 
Benders  cuts  in  the  more  general  framework  of  maximizing  a  submodular  set 
function. 

Branch-and-bound  algorithms  for  the  uncapacitated  plant  location  problem 
use  the  fact  that  it  is  not  necessary  to  constrain  the  variables  y.  .  to  be 

'  J 

integral.  Branching  is  done  on  a  binary  enumeration  tree  with  respect  to  the 
variables  x,.  Bounds  are  obtained  from  one  of  the  two  linear  programming 

J 

relaxations  -  (1.1),  (1.2),  (1.5)  and  (1.3)  or  (1.3').  The  first  algorithms 
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used  the  linear  program  with  (1.3*).  Efroymson  and  Roy  (1966)  showed  that 
this  linear  program  can  be  solved  analytically  so  that  the  bound  at  each  node 
of  the  enumeration  tree  could  be  computed  very  quickly  in  constant  time. 
Improvements  to  Efroymson's  and  Roy's  algorithm  were  made  by  Spielberg 
(1969a),  Khumawala  (1972)  and  Hansen  (1972).  However,  when  (1.3* )  is  used 
in  a  linear  programming  relaxation,  the  bounds  obtained  are  generally  not 
sufficiently  strong  to  curtail  the  enumeration  adequately. 

As  we  observed  previously,  the  constraints  (1.3)  imply  (1.3')  but  not 
conversely.  The  linear  programming  relaxation  that  uses  (1.3)  is  called  the 
string  linear  programming  relaxation  abbreviated  by  SLP.  Revelle  and  Swain 
(1970),  among  others,  observed  that  SLP  is  so  effective  that  its  solution  is 
very  often  integral.  Thus  a  branch -ar.d-bound  algorithm  that  uses  SLP  to 
compute  bounds  is  very  likely  to  perform  well  in  the  sense  that  very  little 
(if  any)  enumeration  will  be  required.  However,  because  of  its  size.  It  is 
not  efficient  to  solve  SLP  directly  by  the  simplex  method. 

Much  of  the  recent  research  on  UL  has  involved  the  development  of 
special  purpose  algorithms  for  solving  SLP.  Marsten  (1972)  used  parametric 
linear  programming  and  a  special  Implementation  of  the  simplex  method. 
Garfinkel,  Neebe  and  Rao  (1974)  used  Dantzig-Wolfe  decomposition.  Schrage 
(1975)  devised  a  generalized  simplex  method  to  treat  the  variable  upper 
bounds  (1.3).  Guignard  and  Spielberg  (1977)  suggested  a  version  of  the 
simplex  method  that  pivots  only  to  Integral  vertices  of  the  polytope  (1.2), 
(1.3),  (1.6).  Cornuejols  and  Thizy  (1982b)  used  a  primal  subgradient 
algorithm. 


uJki 


Dual  algorithms  or  algorithms  that  solve  the  dual  of  the  strong  linear 
programming  relaxation  have  the  advantage  that  upper  bounds  are  obtained  from 
any  dual  feasible  solution.  Thus  a  sufficiently  good  bound  may  be  obtained 
to  fathom  a  node  of  the  enumeration  tree  prior  to  solving  the  dual  to 
optimality. 

Bllde  and  Krarup  (1977)  and  Erlenkotter  (1978)  used  heuristic  methods  to 
obtain  a  near-optimal  solution  of  the  dual.  Erlenkotter  went  a  step  further 
by  using  the  complementarity  slackness  conditions  of  linear  programming  to 
improve  this  bound.  His  procedure  was  so  effective  that  in  45  out  of  the  48 
problems  that  he  tested,  optimality  was  reached  at  the  first  node  of  the 
branch-and-bound  algorithm.  His  OUALOC  code  appears  to  outperform  all  exist¬ 
ing  algorithms. 

A  Lagranglan  dual  of  the  formulation  (1.1)-(1.4),  proposed  by  Geoffrion 
(1974),  is  obtained  by  weighting  the  constraints  (1.2)  by  multipliers  and 
placing  them  in  the  objective  function.  It  can  be  solved  using  subgradient 
optimization,  see  Held,  Wolfe  and  Crowder  (1974).  The  Lagrangian  approach 
can  also  be  used  for  the  p-facillty  location  problem,.  Some  computational 
results  are  reported  in  Narula,  Ogbu  and  Samuelsson  (1977),  in  Cornuejols, 
Fisher,  Nemhauser  (1977b)  and  in  Mulvey  and  Crowder  (1979).  Krarup  and 
Pruzan  (1983)  mention  a  different  Lagrangian  dual  obtained  when  constraints 
(1.3)  (instead  of  (1.2))  are  weighted  by  multipliers  and  placed  in  the 
objective  function. 


3.  Computational  Complexity 

An  algorithm  is  said  to  be  a  polynomial-time  algorithm  for  problem  P, 
if  for  all  instances  of  P  (possible  data  sets),  the  computing  time  of  the 
algorithm  can  be  bounded  by  a  polynomial  function  of  the  data  size.  If  l 
measures  the  data  size  and  k  is  the  order  of  the  polynomial,  we  say  that 
the  computing  time  of  the  algorithm  is  0(L  ).  Sometimes  it  is  more  conve¬ 
nient  to  express  the  computing  time  as  a  function  of  basic  data  parameters, 
such  as  the  dimension  of  a  matrix  or  the  number  of  nodes  in  a  graph.  Then, 
but  only  then,  it  is  assumed  that  all  arithmetric  operations  and  compari¬ 
sons  are  performed  in  unit  time. 

A  fundamental  theoretical  question,  also  of  some  practical  importance, 
is  whether  a  given  combinatorial  optimization  problem  can  be  solved  by  some 
polynomial-time  algorithm.  Oenote  by  P  the  class  of  problems  that  can  be 
solved  in  polynomial -time,  i.e.  by  some  polynomial -time  algorithm.  For  most 
combinatorial  optimization  problems  of  practical  interest,  the  question  - 
are  they  in  P?  -  has  not  been  answered.  A  significant  step  was  made  by  Cook 
(1971)  and  Karp  (1972)  who  introduced  the  notion  of  NP-complete  problems. 

This  is  a  class  of  combinatorial  problems  that  are  equivalent  in  the  sense 
that  either  all  or  none  of  these  problems  can  be  solved  by  a  polynomial-time 
algorithm. 

At  present  no  polynomial-time  algorithm  is  known  for  solving  any 
NP-complete  problem  and  it  has  been  widely  conjectured  that  none  exists.  A 
problem  is  said  to  be  NP-hard  if  the  existence  of  a  polynomial-time  algorithm 
to  solve  it  would  imply  that  all  NP-complete  problems  can  be  solved  by  a 
polynomial-time  algorithm.  Thus  to  show  that  a  problem  (P)  is  NP-hard  it 
suffices  to  find  a  polynomial  transformation  that  reduces  a  known  NP-complete 
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problem,  see  e.g.  the  .comprehensive  list  given  by  Garey  and  Johnson  (1979), 
to  the  problem  (P). 


Theorem  3.1  The  uncapacitated  plant  location  problem  is  NP-hard. 

Proof :  We  need  to  introduce  the  vertex  cover  problem: 

Given  a  graph  G  and  an  integer  k,  find  whether  there  exists  a  subset 
of  k  vertices  of  G  that  cover  all  the  edges  of  G.  (Vertex  v  is  said 
to  cover  edge  e  if  v  is  an  endpoint  of  e.)  The  vertex  cover  problem  is 
NP-complete,  see  Karp  (1972)  or  Garey  and  Johnson  (1979).  We  reduce  it  to 
Ul. 

Consider  a  graph  G  *  (V,E)  with  vertex  set  V  and  edge  set  E. 
Construct  an  instance  of  UL  with  the  set  of  potential  facilities  J  *  V  and 
set  of  clients  I  *  E.  Let  c^  »  0  for  all  i  e  E  and  j  e  V,  and  let 
f.  *  1  for  all  j  e  V.  This  transformation  is  polynomial  in  the  size  of  the 

J 

graph. 

Note  that  the  instance  of  UL  defined  in  this  way  consists  of  covering 
all  the  edges  of  the  graph  G  with  the  minimum  nimber  of  vertices.  Thus  an 
optimal  solution  of  UL  provides  the  answer  to  the  vertex  cover  problem.  This 
proves  that  UL  is  NP-hard.  0 

A  polynomial  transformation  that  reduces  a  known  NP-hard  problem  to  a 
problem  (0)  shows  that  (0)  is  also  NP-hard.  An  immediate  corollary  of 
Theorem  3.1  is  that  the  p-facility  location  problem  is  NP-hard  since  solving 
it  for  every  p  *  l,...,n  provides  a  solution  to  UL. 

Although  UL  is  NP-hard,  some  special  cases  can  be  solved  in  polynomial- 
time.  Kolen  (1982)  has  shown  that  UL  is  solvable  in  time  0(r  )  when  the 
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problem  Is  defined  as  a  tree  with  r  nodes  and  certain  other  assumptions  are 
satisfied.  These  are  that  the  clients  as  well  as  the  facilities  are  located 
at  nodes  of  the  tree.  A  length  Is  associated  with  each  edge  of  the  tree  and 
d^j  Is  the  distance  between  nodes  1  and  j.  Kolen  solves  the  formulation 
(1,1)’,  (1.2)— (1.4)  and  shows  that  the  strong  linear  programming  relaxation 
always  has  an  Integral  optimal  solution. 

Another  Interesting  case  of  the  uncapacitated  plant  location  problem 
that  can  be  solved  In  polynomial  time  was  discovered  by  Krarup  and  Bllde 
(1977).  In  this  Instance  too,  the  crux  Is  that  the  strong  linear  programming 
relaxation  always  has  an  Integral  optimal  solution.  The  conditions  required 
by  Krarup  and  Bilde  generalize  those  obtained  when  a  classical  economic  lot 
size  problem  is  formulated  as  an  uncapacitated  plant  location  problem. 

Finally,  Barany,  Edmonds  and  Wolsey  (1983)  have  given  a  polynomial -time 
algorithm  for  a  tree  partitioning  problem  that  contains  both  Kolen's  and 
Krarup 's  and  Bilde's  problems. 


Suppose  we  are  given  a  feasible  solution  to  UL  that  is  claimed  to  be 
optimal  or  nearly  optimal  (within  a  specified  absolute  or  relative 
tolerance).  We  know  of  only  two  ways  to  verify  this  claim. 

a.  enumeration:  compare,  perhaps  implicitly,  the  value  of  this  feasible 
solution  to  all  others. 

b.  bounding:  determine  an  upper  bound  on  the  optimal  value  of  all  feasible 
solutions  that  is  sharp  enough  to  verify  the  claim. 

Enumeration  is  useful  algorithmically  only  when  it  can  be  done  implic¬ 
itly.  Generally,  this  means  that  the  emmerative  approach,  as  in  a 
branch -and -bound  algorithm,  uses  upper  bounds  to  curtail  the  enuneration. 
Conversely,  an  algorithm  whose  primary  thrust  is  bounding  may  need  to  resort 
to  some  enumeration  to  verify  the  claim. 

The  point  is  that  good  upper  bounds,  as  well  as  good  feasible  solutions, 
are  crucial  in  solving  UL,  as  for  that  matter,  any  hard  combinatorial 
optimization  problem.  We  will  see,  however,  that  UL  has  many  features  that 
make  it  a  relatively  easy  NP-hard  problem. 

Duality  plays  a  key  role  in  the  determination  of  upper  bounds.  The  dual 
of  the  strong  linear  programming  relaxation  given  by  (1.1)-(1.3),  (1.5)  is 

(4.1)  W  *  min  l  u,  +  £ 

id  1  J 

(4.2)  Uj  +  w.j  2.  cij  3,1  i  «  I,  j  €  J 

(4.3)  ->  w^  +  t j  >  -f.  all  j  c  J 


(4.4) 


t,  >  0 


all  i  €  I,  j  £  J. 


We  can  eliminate  variables  and  constraints  from  this  formulation  by 
noting  that: 

a.  for  given  w^,  (4.1)  implies  that  we  would  like  to  make  tj  as  small 

as  possible.  Thus  (4.3)  and  (4.4)  imply  that  t.  •  (  £  w..  -  f, )+, 

j  ‘j  J 

where  (a;  »max(0,a). 

b.  for  given  u^,  (4.1)  also  implies  that  we  would  like  to  make  w^  as 
small  as  possible.  Thus  (4.2)  and  (4.4)  imply  that  w-^  »  (c^  -  u^*. 

Thus 

(4.5)  tj  *  [  l  (c^  -  u^*  -  fj]+  all  j  e  J. 

If  we  think  of  the  u- 's  as  prices  associated  with  the  clients,  the  t^'s 
are  the  profits  from  the  facilities  relative  to  these  prices.  In  other  words, 
if  someone  agreed  to  pay  us  u^  for  the  right  to  serve  the  ith  client, 
we  would  be  willing  to  sell  the  jth  facility  for  the  price  t^  given  by  (4.5). 
Substituting  (4.5)  into  (4.1)  yields  the  condensed  dual 

(4.6)  W  *  min{  £  ui  +  l  [  l  (c^  -  u- )+  -  fj+). 

u  i«I  1  jeJ  id  1J  1  J 

The  dual  is  then  to  determine  a  minimun  sum  set  of  prices  u^  for  the 
clients  and  consequently  a  minimun  sun  set  of  prices  for  the  facilities  so 
that  we  would  agree  to  sell  the  operation.  It  tells  us  the  linear  program¬ 
ming  approximation  of  the  worth  of  our  assets. 

If  I  [(Clj  -  u/  -  fj]+  >  0,  then  some  can  be  increased  without 

increasing  the  objective  function  (4.6).  Also  if  u-  >  max  c^,  then  u, 

1  jeJ  1J  1 


can  be  decreased  without  increasing  the  objective  function  (4.6).  These 
observations  yield  a  second  condensed  dual 


(4.7) 

W  ■  min  >  u^ 
u  i«I  1 

(4.8)  I  (c 

id  iCij 

1  ‘  ui>*  '  'ii° 

all 

j  «  J 

(4.9) 

u,  <  max  c,j 

1  “  jeJ 

all 

i  «  I. 

Although  both  condensed  dual  formulations 

are 

nonlinear 

important  because  they  contain  only  m  variables;  furthermore  for  any  value 
of  u,  we  obtain  the  upper  bound 


(4.10) 


W(u)  »  I  ui  ♦  l  [  l  (cH  -  u  )+  -  fj+. 
iel  1  j«J  id  1  J 


When  (4.8)  holds,  the  upper  bound  reduces  to 

(4.11)  W(u)  *  l  u,. 

id  1 

A  Lagrangian  dual  of  (1.1)-(1.4)  is  obtained  by  weighting  the  constraints 
(1.2)  by  multipliers  u^  and  placing  them  in  the  objective  function.  Let 

(4.12)  l(u)  *  max[  l  l  cHyii  -  l  f^.  +  l  u^l  -  l  y,.)] 

id  jeJ  jcJ  J  J  id  1  j€J  J 

subject  to  (1.3)  and  (1.4). 
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Frequently,  a  Lagrangian  dual  provides  a  tighter  upper  bound  than  a 
linear  programming  dual.  But  this  is  not  the  case  here. 

Proposition  4.2  L(u)  *  W(u)  for  all  u. 

Proof:  L(u)  «  max  £  'i  (c - 4  -  u. )y. .  -  l  f4x4  +  l  u-  subject  to 

id  j«J  1J  1  1J  jc J  J  J  id  1 

(1.3)  and  (1.4).  Hence  y^  *  Xj  if  c^  •  -  u^  >  0,  y^.  »  0  if  c^.  -  u^  <  0 

and  y^  *  0  or  x^.  if  c^  -  *  0.  Thus 


(4.13) 


L(u)  *  max  l  [  I  (Cj.  -  uj*  -  f-]  x,  +  \  u. 

x.<(0,l>  jeJ  id  1J  1  J  J  id  1 


l  t,  +  l  u-  *  W(u) 
jtJ  J  id  1 


where  t.  is  given  by  (4.5).  This  is  true  since  x,  *  1  if  '[  (c, .  -  u- )  + 
J  J  id  'J 

f.  >  0,  x.  =  0  if  l  (c • 4  -  u- )+  -  f ,  <  0  and  x.  *  0  or  1  otherwise.  0 
J  J  ‘J  *  J  J 

Corollary  4.3  L  3  min  L(u)  *  W. 

u 

When  u  satisfies  the  constraints  (4.8),  the  solution  of  (4.13) 
satisfies  the  complementarity  conditions 


(4.14)  (  l  (cH  -  uj+  -  f,)x,  -  0  all  j  e  J. 

iel  J  J  J 


Eauation  (4.14)  suggests  that  if  u  satisfies  (4.8),  to  find  a  good  primal 
solution  we  should  only  consider  ooening  those  facilities  for  which 


5.  Heuristics 


The  combinatorial  formulation  of  the  uncapacltated  plant  location  problem 

max  z(S)  where  z(S)  *  T  max  c4,.  -  7  f4  can  be  viewed  as  a  condensed, 
ScJ  1«I  jeS  S  J 

nonlinear  primal  that  depends  only  on  the  values  of  the  sets  S  £  J.  Based 
on  this  observation  numerous  authors  have  proposed  heuristics  that  Iterate 
on  the  set  S  of  open  locations  and  avoid  an  explicit  Integer  programming 
formulation.  Two  of  the  most  basic  heuristic  approaches  are  described 
below:  the  greedy  and  Interchange  heuristics. 


The  Greedy  Heuristic. 

Start  with  no  facilities  open.  Given  a  set  S  of  open  facilities,  open 
that  facility  j  f  S  whose  Incremental  value  p.(S)  *  z(S  u  {j})  -  z(S)  Is 
as  large  as  possible,  and  Is  positive.  If  no  such  facility  exists  stop  with 
the  set  S  of  open  facilities. 

Formally 


Initialization 
Iteration  t 


S°  ’  '•  P*M  ■  I.1!  '  'J'  '  ■  '• 

Find  j*  ■  arg  max  p^S*"*). 

,  ct-l  J 

jfs 

If  pj  (S*-1)  £  0  and  t  >  I,  stop.  The  set  St_1  is  the 

greedy  solution  with  value  ZS  ■  z(St’^).  If  t  *  1,  the 
greedy  solution  Is  S1  *  {jj}. 

If  pjt(st"1}  >  °»  st  "  st'L  u  ^t}- 


Set  t  «■  t  +  1. 


The  greedy  heuristic  requires  at  most  n  iterations  and  each  Iteration 

2 

requires  0(nm)  calculations.  Thus  the  overall  running  time  Is  0 (n  m). 


We  wil 1  use  the  following  small  example  to  Introduce  and  motivate  the 
Ideas  developed  subsequently.  Real-world  problems  are  typically  much  large 
(e.g.  m  ■  n  »  100).  Consider  the  uncapacitated  facility  location  problem 
defined  by  the  data: 


m  *  4,  n  »  6,  f  *  (3, 2, 2, 2, 3, 3)  and 


Applying  the  greedy  heuristic  to  the  example  yields 

iteration  1:  (Pl(#> . P6(<b) )  »  (16,15,15,12,10.13). 

Hence  j  ^  »  1  and  S1  »  (1). 

Iteration  2:  (p2({1» . p6«l»)  ■  (1, 0,-1, -1,-1). 

Hence  j2  *  2  and  S2  -  {1,2}. 

Iteration  3:  (p2({l»2}),.*.,pg({l,2}))  ■  (0,-1, -1,-1). 

?  c  p 

The  set  S  of  value  Z  *  z( S  )  ■  17  is  the  greedy  solution. 

We  can  now  use  (4.10)  to  obtain  upper  bounds  on  the  values  of  the  feasi 
ible  solutions  produced  by  the  greedy  heuristic.  In  fact,  such  an  upper 
bound  can  be  associated  with  any  S  c  J. 


Define 


Note  that 


TT.(S)  *  max  c.  .  all  1  «  S. 
jeS  1J 


z(s)  -  l  «Tt(s)  -  i  f. 

i«i  1  us  3 


Thus,  from  (4,10) 


since 


P.(S)  *  I  (c. .  -  IT.  (S))+  -  f .  all  j  f  S. 

J  U  i  J 


W(u(S))  -  l  u\(S)  *  l  [Pi(S)]+ 
1  € I  1  US  3 


c^.  -  u^(S)  <^0  for  all  1  «  I  and  j  €  S. 


In  particular  If  S  Is  the  final  set  chosen  by  the  greedy  heuristic 

then,  by  the  stopping  criterion,  S ^ (SG)  £0  j  ^  SG  so  that  wCu(SG))  * 

"  u^  (S  ).  We  have  shown  that  the  greedy  solution  deviates  from 

optimality  by  at  most  7  f.,  which  suggests  that  It  will  yield  a  small 

j«sG  3 

error  when  the  fixed  costs  are  small  in  comparison  to  the  profits. 
Furthermore,  we  may  obtain  a  better  bound  by  considering  all  of  the  sets 

produced  by  the  greedy  algorithm.  Let  u?  *  min  c^  all  i  «  I  and 

j«0  13 

u*  ■  u(Sk),  k  »  l,...,t-l.  Define  a  dual  greedy  value  by  WG  *  min  W(u^). 

k 

In  the  example,  17^  *(0  0  0  0),  W (TT^ )  ■  "  p  (<>)  *  83,  IT1  »  (6, 6, 5, 2), 

JeJ  3 


WOO  *  z  lit*  1  -  20,  u‘  >  (6, 8,5, 3)  and  W(uc)  »  7  u7  *  22.  Hence 
i«I  1«I  1 

WG  -  20. 

The  bound  we  have  given  so  far  is  an  a  posteriori  bound  for  a  particular 
instance  of  UL.  In  fact,  a  general  relationship  between  ZG  and  WG  that 
a  priori  applies  to  all  Instances  of  UL  is  given  by 

Theorem  5.1  [Cornuejols,  Fisher  and  Nemhauser  (1977b)] 

Za  >  (®ji)WG  *  (I)R 


where  e  is  the  base  of  the  natural  logarithm  and 


'  m1n  cii  “  I  V 

jcJ  J 


Ul  je J 


A  proof  of  Theorem  5.1  that  uses  linear  programming  duality  is  given  in 
Fisher,  Nemhauser  and  Wolsey  [1978], 

For  the  p-facility  location  problem  with  c.  .  >_  0  for  all  i  and  j 

1  J 


and  f.  *  0  all  j  e  J,  we  have  R  j>  0.  Thus  we  achieve  a  simple  data 


Independent  statement  of  Theorem  5.1. 


Corollary  5.2  [Cornuejols,  Fisher  and  Nemhauser  (1977)]  For  the  p-facility 
location  problem  with  c,  .  _>  0  for  all  1  and  j  and  f  .  »  0  all  j  £  J 

*  J  J 


There  are  families  of  p-facility  location  problems  for  which  this  bound 
is  achieved  asymptotically.  Furthermore,  since  ZG  <  Z  <  W  <  WG 


•**{ r  •  I  •  fV)  •  °-86- 
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There  are  several  variations  and  generali zations  of  the  greedy  heuristic 
for  which  bounds  similar  to  those  of  Theorem  5.1  and  Corollary  5.2  are  known, 
see  Cornuejols,  Fisher  and  Nemhauser  [1977b]  and  Nemhauser,  Wolsey  and  Fisher 
[1978].  For  example,  we  can  begin  with  all  facilities  open  and  at  each 
iteration  close  a  facility  that  gives  the  largest  improvement  in  the  objec¬ 
tive  function  so  long  as  such  a  facility  exists.  A  generalization  of  the 
greedy  heuristic  is  to  start  with  the  family  consisting  of  all  sets  of  cardi¬ 
nality  k,  for  some  fixed  k,  and  apply  greedy  to  each  of  these  (£)  initial 
sets  separately;  we  then  choose  the  best  of  the  resulting  (")  solutions. 

None  of  these  variations  or  generalizations,  however,  improve  the 
worst-case  bound  of  the  greedy  heuristic.  In  fact,  what  is  remarkable  about 
the  bound  on  the  greedy  heuristic  Is  not  its  value,  but  that  no  polynomial¬ 
time  procedure  of  any  degree  whatsoever  is  known  for  p-facility  location 
problems  that  has  a  better  worst-case  performance. 

The  salient  feature  of  the  greedy  heuristic  is  that  the  maximum  possible 
improvement  is  made  at  each  step.  If  this  is  not  done,  worst-case  perfor¬ 
mance  deteriorates,  even  if  a  broader  choice  of  improvements  is  considered. 

An  example  of  such  a  heuristic  is  generalized  interchange,  see  Nemhauser, 
Moisey  and  Fisher  (1978).  Here  we  begin  with  an  arbitrary  set  S°.  Given 
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St-1,  at  Iteration  t  Me  select  any  set  St  such  that  z(St)  >  z(St_1), 

j St\St_1  j  <  1  and  JSt"Kst|  1  or  stop  If  no  such  S*  exists.  Thus, 

at  each  Iteration,  we  are  allowed  to  open  a  facility,  close  a  facility  or  do 

both  so  long  as  an  Improvement  Is  made. 

The  worst«case  bound  of  generalized  Interchange  Is  weaker  than  that  of 

greedy.  For  example,  under  the  conditions  of  Corollary  5.2,  this  heuristic 

can  guarantee  only  to  find  a  solution  of  value  at  least  half  of  the  optimum 

value.  Of  course,  by  starting  with  a  greedy  solution,  the  bounds  of  Theorem 

5.1  and  Corollary  5.2  are  obtained.  Nevertheless,  they  are  not  strengthened 

by  applying  generalized  Interchange.  On  the  other  hand,  greedy  followed  by 

generalized  Interchange  seems  to  give  good  empirical  performance,  see  Hansen 

(1972)  and  Comuejols,  Fisher  and  Nemhauser  (1977b). 

The  dual  solution  7{S)  given  by  (5.1)  is  motivated  by  its  use  in 

obtaining  a  bound  on  the  given  primal  solution  S.  Conversely,  given  r  dual 

solution  that  satisfies  (4.8),  the  complementarity  conditions  (4.14)  suggest 

considering  a  primal  solution  In  which  x.  *  0  if  T  (c. .  -  u J+  -  f .  <  0 

J  1j  i  1 

Let 

J(u)  -  (j:  :  (Cii  -  ui )+  -  f .  *  0). 

lei  1  J 

The  best  solution  that  satisfies  complementarity  is  obtained  by  solving 

max  {  "  max  c4 ,  -  "  f,}, 

ScJ(u)  1«I  jeS  jeS  J 


but  this  problem  may  not  be  much  easier  to  solve  than  UL  itself.  Instead,  we 
take  any  minimal  set  K(u)  c  J(u)  that  satisfies 


(5.2) 


max  c-,.  *  max  c. .  all  i  e  I. 
jeK(u)  13  JcJ(u)  13 


Proposition  5.3  Given  a  u  that  satisfies  (4.8)  and  u^  <  max  c^.. 

-  '  “  J«J(u>  'J 

all  i  e  I,  and  a  primal  solution  K(u)  defined  by  (5.2),  let 

ki  s  |{j  e  K(u):  >  u -  > | .  If  k-  £  1  all  i  e  I,  then  u  is  an  optimal 

set  of  open  facilities. 


Proof: 


If  ki  *  0 


z(K(u))  *  l  max  c. .  -  l  fv 

iel  jeK(u)  13  jeK(u)  3 


max  c,,  *  u,  *  u,  +  >  (c-.  -  u-)+ 

JeK(u)  13  1  1  13  1 


jeK(u) 


and  if  k^  3  1 


max  c. .  =  u,  +  l  (ci1  -  uj+. 
JeK(u)  13  1  jcK(u)  13  1 


Hence,  if  k^  £  1  all  iel. 


z(K(u) )  *  l  I  (cH  -  u.)  -  l  f.  +  l  u, 

iel  jeK(u)  13  1  jeK(u)  3  iel  1 


jeK(u) 


‘l  (  L  (cn  -  ui)+  -  f.)  +  l  ui 
jeK(u)  iel  3  3  iel 


•l  ui  *  W(u)  by  (4.11).  3 

iel  1 


In  our  example,  with  u  ■  (6  6  4  3),  we  obtain  J(u)  ■  K(u)  ■  {2,3,4} 

and  (k^.kj.kj.k^)  *  (1  1  I  0).  Hence  these  are  optimal  primal  and  dual 

solutions  of  value  T  u4  »  19. 

l2l  1 

While  Proposition  5.3  may  permit  us  to  recognize  an  optimal  solution,  it 

is  limited  to  those  cases  In  which  min  W(u)  ■  Z  and  even  then  It  Is  still 

u 

necessary  to  find  an  appropriate  u  and  K(u). 

Dual  Descent 

Dual  descent  Is  a  heuristic  that  begins  with  a  u  satisfying  (4.8)  and 
then  attempts  to  decrease  the  u^'s  one-at-a-time  while  maintaining  (4.8), 
see  Erlenkotter  (1978)  and  Bilde  and  Krarup  (1977).  It  is  surprisingly 
effective,  but  not  fail  safe.  In  finding  a  u  that  satisfies  the  conditions 
of  Proposition  5.3.  This  descent  approach,  with  some  embellishments,  is  the 
inner  loop  of  Erlenkotter's  0UAL0C  algorithm.  The  basic  descent  method 
precedes  as  follows: 

Begin  with  u?  *  max  c. .  all  i  e  I.  Cycle  through  the  indices  i  c  I 
1  j£J 

one-by-one  attempting  to  decrease  to  the  next  smaller  value  of  c^j. 

If  one  of  the  constraints 

(5.3)  2  <cii-ui)+<fi  311 

1«I  ”  J 

blocks  the  decrease  of  u^  to  the  next  smaller  c^,  is  decreased  to  the 
minimum  value  allowed  by  the  constraint.  When  all  of  the  u^'s  are  blocked 
from  further  decreases,  the  procedure  terminates. 

The  reason  for  decreasing  only  to  the  next  smaller  c^  . ,  rather  than 
to  the  smallest  permissible  value,  is  to  keep  the  k.  of  Proposition  5.3  as 
small  as  possible. 


Applying  dual  descent  to  our  example  yields  the  results  shown  in 
Table  5.1.  For  the  first  four  steps,  each 


Table  5.1 


of  the  u^s  is  decreased  to  the  second  max  in  the  row.  Now  is  con¬ 
sidered  again,  but  cannot  be  decreased  because  the  constraints  (5.3)  for 
j  =  3  would  be  violated.  Similarly,  a  decrease  of  u^  would  violate  (5.3) 
for  j  *  2.  However,  u^  can  be  decreased;  but  it  is  decreased  only  to  4 
because  (5.3)  becomes  active  for  j  *  4  when  u^  =  4.  Finally  cannot 
be  decreased  because  of  (5.3)  for  j  *  2.  This  completes  the  dual  descent 

with  u  «  (6  6  4  3)  and  W(u)  *  7  u4  *  19.  Now,  as  noted  above,  J(u)  * 

iel  1 

K(u)  *  (2,3,4)  and  we  also  obtain  a  primal  solution  of  value  19. 

A  possible  improvement  of  dual  descent,  which  is  likelier  to  produce  a 
primal  and  dual  pair  for  which  Proposition  5.3  applies,  is  obtained  by 
modifying  the  order  in  which  the  u^'s  are  considered  as  candidates  to 
decrease.  In  particular,  rather  than  just  cycling  through  the  u-'s,  let 
Q^u)  *  (j:  c^  -  Uj  >_  0).  Then  if  is  decreased  and  descent  terminates 


k .  £  jQ-(u) j .  Hence  we  choose  u$  next  if  |0$(u)|  £  |Q^(u) |  for  all 

i  «  I. 

Suppose  dual  descent  terminates  with  the  dual  solution  u*  and  we 
determine  a  primal  solution  given  by  K(u*)  such  that  Proposition  5.3  fails 
to  verify  optimality.  Then  there  exists  an  i  such  that  ki  >  1 .  By 
increasing  u^  to  its  previous  value  and  reapplying  dual  descent,  it  may  be 
possible  to  improve  the  dual  solution  further. 

We  have  sketched  all  of  the  basic  steps  of  Erlenkotter's  heuristic.  An 
example  in  which  it  does  not  find  an  optimal  solution  to  (4.6)  is  given  by 

/2  2  0  \ 

(5.4)  f  *  (2  2  2),  C  *  f  2  0  2  j  . 

\0  2  2/ 

In  this  instance  of  UL,  beginning  with  u®  *  (2  2  2),  dual  descent  terminates 
with  u  *  (0  2  2)  and  the  embellishments  don't  help.  However  an  optimal  u 
is  (1  1  1)  yielding  W  »  3  and  1*2.  Nevertheless,  this  heuristic  has 
performed  extremely  well  on  the  problems  that  Erlenkotter  has  considered. 

He  reports  that  in  45  of  48  problems  tested,  the  heuristic  found  an  optimal 
solution.  To  provide  the  capability  of  finding  an  optimal  solution  and 
proving  optimality,  the  heuristic  is  imbedded  in  a  branch-and-bound  algorithm 
called  DUALOC.  Given  its  simplicity,  speed  and  availability,  DUALOC  may  be 
the  most  efficient  way  to  solve  Ul.  However,  it  could  bog  down  on  hard 
problems  in  which  the  heuristic  bound  is  not  as  good  as  the  linear  program¬ 
ming  bound.  Thus  one  is  motivated  to  develop  efficient  algorithms  for 
solving  the  strong  linear  programming  relaxation. 
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6.  Algorithms  and  Reformulations  of  the  Strong  Linear  Programming  Relaxation 
The  strong  linear  programming  relaxation  (SLP)  of  the  uncapacitated 
facility  location  problem  (UL)  is 


(6.1) 

ZL?  «  max  1 
i  el 

n 

J«J 

Cijv '  jL  Vi 

(6.2) 

n 

j*J 

yij 

*  i 

all 

i  e  I 

(6.3) 

yij 

-  x ,  <_  0 

all 

i  cl,  j  £  J 

(6.4) 

yij 

>  o,  Xj  >  0 

all 

i  e  I,  j  e  J. 

For  general  integer  programs,  a  linear  programming  relaxation  must  be 
used  in  conjunction  with  cutting  planes  or  enuneration  to  obtain  an  optimal 
integral  solution.  However,  for  reasons  which  are  barely  understood,  SLP  is 
an  unusually  powerful  relaxation  of  UL  in  the  following  sense. 

Observation  6.1  Very  frequently,  SLP  has  an  optimal  integral  solution. 

Of  course,  it  is  not  true  that  all  of  the  extreme  points  of  the 
polyhedron  (6.2)-(6.4)  are  integral.  Example  (5.4)  has  a  unique  optimal 
solution  with  x  *  (1/2  1/2  1/2)  and  it  is  easy  to  construct  infinite 
families  of  objective  functions  for  which  the  unique  optimal  solution  to 
(6.1)-(6.4)  is  fractional.  Nevertheless,  randomly  constructed  objective 
functions  and  the  few  encountered  in  practice  that  have  appeared  in  the 
literature  strongly  support  Observation  6.1. 


A  challenging  problem  for  the  combinatorial  mathematician  is  to  make 
Observation  6.1  precise.  For  the  practitioner.  Observation  6.1  means  that 


an  efficient  method  for  solving  SLP  win  also  be  an  efficient  method  for 
solving  most  Instances  of  UL.  Even  If  one  Is  unlucky,  a  very  good  upper 
bound  Is  typically  attained  so  that  a  branch-and-bound  algorithm  should 
terminate  rapidly. 

Even  relatively  modestly  sized  Instances  of  UL  cannot  be  solved  by  a 
standard  mixed  Integer  programming  package  that  uses  SLP  as  a  linear 
programming  relaxation  because  of  the  large  number  of  constraints  (6.3).  For 
example,  m  *  n  *  100  yields  a  problem  with  more  than  10,000  constraints.  In 
this  section,  we  consider  several  approaches  for  solving  very  large 
structured  linear  programs.  We  will  begin  by  briefly  mentioning  two  direct 
approaches  to  eliminating  the  difficulty  caused  by  the  large  number  of 
constraints  (6.3).  Then  we  will  apply  some  well-known  reformulations  and 
algorithms  including  Lagrangian  duality  and  subgradient  optimization, 
Dantzlg-Wolfe  decomposition,  8enders  decomposition  and  subgradient 
optimization  on  the  primal.  Connections  among  those  approaches  will  be 
noted.  Finally,  we  will  consider  a  reformulation  that  involves  the  reduction 
of  the  matrix  C  into  an  interesting  canonical  form. 


Direct  Approaches 

The  constraints  y, .  _<  x.  are  generalizations  of  simple  upper  bound 

•  J  J 

constraints  in  which  the  upper  bounds  themselves  are  variables.  It  is 
well-known  how  to  handle  fixed  upper  bound  constraints  in  the  simplex  method 
without  expanding  the  dimension  of  the  basis  to  include  them.  Schrage  (1975) 
has  generalized  this  idea  to  incorporate  variable  upper  bounds.  He  reports 
computational  results  obtained  by  applying  his  method  to  SLP. 


An  alternative  is  to  generate  the  constraints  y. .  £  x.  as  cuts  only 
when  they  are  violated  in  an  optimal  solution  to  the  weak  linear  programming 
relaxation.  This  Idea  has  been  tested  by  Morris  (1978).  In  the  example  of 
the  last  section,  an  optimal  solution  to  the  weak  linear  programming  relax¬ 
ation  Is  %2  *  x3  ■  x4  *  x5  *  1/4  and  y^3  *  y22  *  y34  *  y45  *  1.  We  could 
then  add  the  four  violated  variable  upper  bound  constraints  and  continue  to 
solve  the  linear  program. 

The  direct  approaches  are  primal.  Dual  methods  may  be  superior  for  two 
reasons.  First,  if  SLP  Is  Incorporated  in  a  branch-and-bound  algorithm,  it 
may  not  be  necessary  to  solve  SLP  to  optimality  at  every  node  of  the 
enumeration  tree,  i.e.  at  some  of  the  nodes,  a  dual  feasible  solution  may 
suffice  to  bound  the  subproblem.  Second,  as  we  have  already  shown  in  the 
last  section,  we  can  easily  generate  an  integral  primal  solution  from  each 
dual  solution. 

Lagrangian  Duality  and  Subgradient  Optimization 

The  Lagrangian  L(u)  of  (4.12)  forms  the  basis  of  a  subgradient 

algorithm  for  solving  the  dual  of  SLP.  The  subgradient  algorithm  solves  the 

problem  min  L(u). 
u 

The  function  L(u)  given  by  (4.12)  is  the  maximum  of  a  finite  number 
of  linear  functions.  Therefore  L(u)  is  piecewise  linear  and  convex.  Sub- 
gradient  optimization  [Held,  Wolfe  and  Crowder  (1974)]  has  proved  to  be  a 
useful  method  for  minimizing  unconstrained  piecewise  linear  convex  functions. 
This  approach  is  an  extension  of  the  gradient  method  for  minimizing  smooth, 
nonlinear  convex  functions.  Since  gradients  do  not  exist  at  non-di fferentiabl 


32 


points  of  L(u),  the  gradient  direction  Is  replaced  by  a  subgradient  direc¬ 
tion,  which  will  be  explained  below. 

Given  u^,  an  Iteration  of  the  subgradient  algorithm  generates  a  new  dual 
solution  by  the  formula 

(6.5)  ut+1  ■  u 1  -  ytaL(ut) 

where  dL(ut)  is  a  subgradient  at  u*  and  yt  is  the  stepsize.  If 
3L(ut)  *  0,  then  ut  is  an  optimal  dual  solution. 

Suppose 


L(u) 


r  i 

lei  jeJ 


y 

jcJ 


f  4X 


n 


+  V 

ijl 


where  {x*,  y.n  are  defined  in  the  proof  of  Proposition  4.2.  If  the  (y . r> 

J  1 J  1  J 

are  unique  then 


(6.6) 


altuL  *  1  -  l  ytf  all  1  e  I 
1  jeJ  1J 


Is  the  gradient  of  l(u)  at  u.  However  if  the  }  are  not  unique, 

then  any  direction  given  by  (6.6)  or  convex  combinations  of  such  directions 
is  a  subgradient  direction.  Although  a  step  in  a  subgradient  direction  does 
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£ 


not  guarantee  a  decrease  in  l(u),  it  can  be  proved  that  with  an  appropriate 
choice  of  stepsize  the  iterates  given  by  (6.5)  converge  to  an  optimal  solution 
[Polyak  (1969)].  Cornuejols,  Fisher  and  Nemhauser  (1977b)  have  solved  the 
Lagrangian  dual  by  subgradient  optimization  and  report  computational  results. 

In  our  example,  if  we  start  with  u°  «  (8  8  6  4),  then  L(u°)  *  26  and 
aL(u)  »  (1  1  1  1).  With  a  stepsize  of  y°  *  2,  we  obtain  u1  »  (6  6  4  2) 
and  L(u^)  *  19.  The  solution  to  (4.13)  is  not  unique;  however,  the  solution 
x2  -  x4  «  1,  x.  «  0  otherwise  and  y12  *  y22  «  y34  *  y42  »  1,  yjj  »  0 

otherwise  yields  aL(u*)  *  (0  0  0  0)  and  verifies  the  optimality  of  u1. 


Dantzig-Wolfe  Decomposition 

D 

For  all  non-empty  R  £  I,  let  x.  *  1  if  facility  j  serves  only  these 

J 

clients  In  the  set  R  and  Xj  *  0  otherwise.  If  *  1,  facility  j  yields 

a  profit  of  T  c..  -  f..  Thus  UL  can  be  reformulated  as  the  Integer  program 
ieR  J 


(6.7) 


Z  *  "*x  l  Id  c„  -  f.)x" 

j«J  Rcl  1€R  ^  J  J 


(6.8) 


(6.9) 


?  r  X,  *  1 
jcd  Rii  J 


"  <  1 
Rcl  J  “ 


all  1  £  I 


all  j  £  J 


(6.10) 


£  {0,1}  all  R  C  I,  j  £  J. 

J 


The  equations  (6.8)  state  that  each  client  is  served  by  exactly  one  facility 
and  the  inequalities  (6.9)  state  that  each  facility  can  serve  only  one  set  of 
clients. 


We  are  going  to  study  the  linear  programming  relaxation  of  this  integer 
program  where  (6.10)  is  replaced  by 

(6.11)  _>  0  all  Rc  I,  j  e  J. 

Proposition  6.2  Let  7  be  the  optimal  value  of  the  linear  program  (6.7), 
(6.8),  (6.9)  and  (6.11).  Then  7  *  Z^. 

Proof  Apply  Dantzig-Wolfe  decomposition  to  the  linear  program  (6.1)-(6.4) 
with  master  constraints  (6.2)  and  subproblem  constraints  (6.3),  (6.4)  and 
Xj  <  1  all  j  c  J.  Substituting  the  subproblem  extreme  points  into  (6.2) 
yields  (6.8),  and  (6.9)  are  the  convexity  constraints  for  the  subproblems.  0 

Before  considering  an  algorithm,  we  will  make  some  simplifications.  If 

r  %  a 

l  (c- .  -  f  .)  >  i  (c-.  -  f.  ),  then  =  0  in  every  optimal  solution. 
ieR  J  J  ieR  K  * 

Hence  for  each  R,  we  need  only  one  variable,  say  with  price 


dR 


fim  s  max(  2 

jcJ  ieR 


’ij 


Furthermore  the  constraints  (6.9)  are  superfluous.  This  is  true  because  if 
R  i  R'  *  $,  then  (6.8)  guarantees  that  \R  +  _<  1  and  if  R  n  R'  *  a  and 


j(R)  *  j(R')  »  k,  then  >  0  Implies  dR(jRl  >  dR  ♦  dR1. 

Thus  we  can  restate  the  Integer  program  as 

(6.12)  Z  -  max  T  d  x„ 

Rc  I  K  R 

(6.13)  7  XD  ■  1  all  1  e  I 

R»i  K 

(6.14)  xR  £  {0,1}  all  R  =  I. 

This  formulation  has  2ffl“*  variables  and  m  constraints.  Curiously,  n 
does  not  enter  into  the  size  of  the  problem  so  that  for  fixed  m,  this 
formulation  can  be  solved  In  polynomial -time. 

The  linear  program  (6.12),  (6.13)  and 

(6.15)  xR  >  0  all  R  c  i 

can  be  solved  by  column  generation.  Suppose  we  begin  with  m  columns,  say 
those  for  R  =  {i},  i  *  l,...,m.  Then  the  primal  solution  is  \R  *  1  for 

R  *  {i},  i  *  l,...,m,  and  the  dual  solution  is  u?  *  max(c. .  -  f.)  all 

1  jeJ  1J  J 

i  e  I. 

We  now  see  if  any  of  the  nonbasic  columns  have  a  positive  price.  This 
can  be  done  at  iteration  p  by  solving  for  each  j  £  J  the  subproblem 


(6.16) 


tP  *  max  "  (c. .  -  uf )y.  .  -  f -x  . 


We  obtain 


since  £  (c,.  -  u?)  -  is  the  price  of  variable  x1?. 
i«R  J  J  J 

t?  a  (  £  (c. .  -  u?)+  -  f. )+.  If  t?  »  0  all  j  {  J  then  the  current 

J  'J  *  J  J 

solution  is  optimal.  For  each  k  such  that  tp  >0,  let  be  any  set 
satisfying  {i:  c^k  -  up  >  0}  _c  R^  c.  (i :  c-k  -  up  >_  0} .  We  now  add  to  the 
linear  program  the  variables  xR  for  all  k  such  that  t£  >  0. 

Garfinkel,  Neebe  and  Rao  (1974)  have  obtained  computational  experience  with 
this  type  of  algorithm. 

An  important  feature  of  this  approach  is  that  both  lower  and  upper  bounds 

on  ZLp  are  obtained  at  each  iteration.  By  primal  feasibility,  ZLp  _>  £ 

id 

and  from  (4.12)  we  see  that  Zy,  £  £  t?  +  £  up  *  W(uP).  Moreover,  as 

well  as  obviously  being  a  primal  method,  it  is  also  a  dual  method  that  can  be 
compared  with  solving  the  Lagrangian  dual  by  subgradient  optimization.  Here 
the  u^'s  at  each  iteration  are  determined  by  solving  a  linear  program,  while 
in  the  previous  method  the  dual  variables  are  determined  by  moving  in  the 
direction  of  a  subgradient  of  the  function  l(u).  An  advantage  of  the  colunn 
generation  approach  is  that  its  lower  bounds  are  determined  in  a  less  ad  hoc 
fashion.  Furthermore,  whenever  the  linear  program  has  an  integral  solution  a 
feasible  solution  to  the  integer  program  is  also  found. 

In  our  example,  we  start  with  an  initial  basis  consisting  of  the  four  unit 
columns  R^  *  { i } ,  i  *  1,...,4  with  objective  coefficients  dR  *  dR  *  6, 


X 


d0  a  4  and  dB  a  1.  This  yields  the  dual  solution  u  =  (6  6  4  1). 
k3  k4 

Solving  the  subproblems,  we  obtain  t^  »  (0  2  0  1  0  0)  and  thus 

17  £  ZLp  £  20.  We  generate  two  new  columns  Rg  «  {1  2  4}  and  Rg  =  {3  4} 

with  objective  coefficients  dg  3  15  and  dg  a  6.  The  next  master  linear 

program  yields  the  primal  solution  a  xi  *  1,  X.J  =0  otherwise  and  the 

*5  K3  Ki 

dual  solution  u^  *  (6  6  4  3).  Now  t^  s  0  all  j  e  J  so  that  the  primal 
has  been  solved.  In  terms  of  the  original  variables,  we  have  x^  a  x^  a  1, 
x.  a  0  otherwise. 

J 


Primal  Subgradient  Algorithm 

With  fixed  x.'s,  0  <  x.  <  1  all  j  e  J,  >  x.  >  1,  SLP  is 
J  J  “  jeJ  J  " 


(6.19) 


Z,„(x)  a  -  I  f.x.  +  max  '[  j>  c Hy,, 
LP  jeJ  J  J  i«I  jcJ  1J  1J 


(6.20) 


l  =1  all  i  e  I 


(6.21) 


yij  ^  XJ 


all  wl  j  e  J 


(6.22) 


y ij  >  0  a11  1  £  r*  J  £ 


Let  V-(x)  be  the  profit  from  the  ith  client.  Then  (6 -19) -(6.22)  decomposes 


as  follows 


(6.23) 


ZLp(x)  3  -  l  fjxj  +  l  Vi (x) 
LK  jeJ  J  J  id  ' 


where  for  each  i  e  I 


(6.24) 
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Vi (x)  ’max  jL  CUyij 
(6*25)  l  yu  -  1 

(6.26)  yij-xj  a11  jeJ 

(6.27)  ^ijl0  311  i  «  J* 


This  linear  program  can  be  solved  greedily.  Let  c  .  >  c  •  >  ... 

~  i J 2  ~ 

cu„-  Then  k  *  1—>p  - l*  %  ” 


'ij  3  1  *  l  xi  »  ^ii  55  0 
Jp  k*l  Jk  1Jk 


P-1  p 

otherwise  where  £  x.  <  1  <  £  x4  . 

k-1  Jk  ”  k*l  Jk 

The  dual  of  (6.24)-(6.27)  is  for  each  i  e  I 


> 


(6.28)  V.  (x)  *  min  u.  +  £  x.w.  . 

1  1  jcJ  J  1J 


(6.29)  u.  *  w(j 

^Ci, 

3  •" 

j  C  J 

(6.30)  WjJ 

l  0 

all 

j  £  J. 

To  solve  ( 6 .28) -( 6 .30) ,  observe  that 

"id 

3  (cij 

-  u^)  +  and  that  it 

sufficies  to  consider  e  {c^, ...,cin>.  Hence 


V.(x)  »  min(cik  +  £  x  . ( c , -  c.J+)  all  i  e  I. 
1  kcJ  jeJ  J  1J  1K 


(6.31) 


Since  V^(x)  is  the  minimum  of  a  finite  number  of  linear  functions,  it 
is  a  piecewise  linear  and  concave  function.  Therefore  Z^p(x)  is  also 
piecewise  linear  and  concave  and 


ZLp  -  max  ZLp(x) 


(6.32) 


0  _<  x.  <  1  all  jeJ 

J 


Xj  >  1 


can  be  solved  by  subgradient  optimization.  If  there  Is  a  j  such  that 

c. .  -  f.  >  0,  the  constraint  "  x.  >  1  is  superfluous.  If  not,  we 
icl  13  3  jcJ  3  “ 

can  add  a  constant  M  >  min(f .  -  "  c,.)  to  some  row  of  matrix  C  without 

jcJ  3  icl  13 

changing  the  solution  while  assuring  that  "  x.  >  1  will  be  satisfied  by 

JeJ  3  “ 

any  optimal  solution.  In  the  remainder  of  the  chapter,  we  assume  that 

x.  >  1  is  not  needed. 

3 

A  subgradient  of  Z.  D(x)  at  x  is  of  the  form 


(6.33) 


av(x) 


*  T  (c,  .  -  c.  )+  -  f .  al 
icl  13  ip  3 


where  p  is  determined  from  a  solution  to  the  ith  subproblem  (6.24)-(6.27) , 

i.e.  c.  *  min  (c,  . :  y. .  >  0>.  Note  that  because  of  the  bounds  on  the  vari- 
ip  j  iJ  iJ 

ables,  if  x^  *  0  then  dV(x)j  is  replaced  by  [  ~  (c^.  -  c^p)+  -  fj]+, 
and  if  x.  ■  1  then  aV(x)..  is  replaced  by  min(0,  "  (c<  .  -  c^)4-  -  f.). 


y 
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Cornuejols  and  Thizy  (1982b)  report  their  computational  experience  in 
solving  (6.32)  by  a  subgradient  algorithm. 

The  results  of  attempting  to  solve  (6.32)  by  subgradient  optimization  for 
our  example  are  summarized  in  Table  6.1.  An  optimal  solution  to  SLP  is 
obtained  at  iteration  3,  but  the  subgradient  formula  (6.33)  Is  not  adequate  to 
verify  it. 


iter¬ 

ation 


V  0  11110 

1  0  1/2  1/2  1/2  1/4  0 

0  3/4  1/2  3/4  0  0 

0  11/2  1  0  0 


V(x)  ZLp(x) 


8  8  6  4  17 

7  7  9/2  3  18 

7  15/2  21/4  11/4  37/2 

7  8  6  3  19 


av(x) 


-1  1  0  1-1-1  1/4 
-1  1  0-1  -1  -1  1/4 
-3  -2  0  -2  -2  -2 


This  is  precisely  the  linear  program  that  arises  when  applying  Senders 
decomposition  to  SIP.  But  now  we  have  mn  constraints  of  the  form  (6.35). 


However,  we  can  think  of  these  as  cutting  planes  and  generate  them  only  as  we 
need  them. 

In  particular,  suppose  we  have  only  a  proper  subset  of  the  constraints 
(6.35).  We  solve  the  relaxed  linear  program  and  determine  an  optimal 
solution  (xq,Vq).  Now  we  use  xq  in  (6.24)-(6 .27)  to  determine  V. (xq) 
all  i  €  I.  If 


(6.37) 


Vi(xq)  <  Vq  all  i  £  I, 


then  (xq,Vq)  satisfies  all  of  the  constraints  (6.35)  and  xq  is  an  optimal 
solution.  If  not,  each  i  for  which  (6.37)  is  violated  specifies  a  violated 
constraint  of  the  form  (6.35).  These  are  added  to  the  linear  program  and  we 
continue. 

In  our  example,  we  begin  with  the  constraints  (6.35)  determined  by  the 
second  maximum  in  each  row,  i.e.,  -  2x^  _<  6,  V£  -  2x^  <_  6,  -  x4  £  5, 

and  V4  <  4.  A  solution  to  the  linear  program  is  V*  3  (6  8  5  4)  and  x*  * 
(0  10000).  Then  by  solving  ( 6 .24) -(6.27) ,  we  obtain  V(x^)  =(6803) 
and  generate  the  constraints 


Now  we  obtain  the  solution  V 

.2 


(6863)  and  x  *  (0  1  0  1  0  0).  Since 

2 


V(x  )  *  (6  8  6  3),  all  of  the  constraints  of  (6.35)  are  satisfied  and  x 
is  an  optimal  solution. 

Magnanti  and  Wong  (1977)  have  used  a  variation  of  this  approach  and  have 
developed  stronger  inequalities  in  an  attempt  to  impose  integrality  on  x. 


Canonical  Reduction 


We  now  develop  another  formulation  that  involves  the  disaggregation  and 
aggregation  of  clients,  see  Cornuejols,  Nemhauser  and  Wolsey  (1980).  The 
aggregation  of  two  clients  i^  and  12  means  to  replace  clients  ij  and 
i2  by  a  single  client  1  such  that 


(6.38) 


:1j  *  Ci  jj 


+  c.  .  all  j  £  J  . 
12J 


The  disaggregation  of  client  1  into  two  clients  i^  and  means  to 
replace  1  by  two  clients  i^  and  i2  such  that  (6.38)  holds.  While 
aggregation  is  uniquely  defined,  disaggregation  is  not. 


•.rVrv’ 


yri 


In  general,  aggregation  yields  an  underestimation  of  profit  and  disaggre¬ 
gation  overestimates.  This  is  true,  because  if  (6.38)  is  satisfied  the  greedy 

solution  to  the  linear  program  (6.24) -(6 .27)  implies  V-(x)  <  V.  (x)  +  V-  (x) 

’l  ’2 

for  all  x  such  that  0  £  x^  £  1.  We  say  that  aggregation  or  disaggregation 
is  val id  when  for  all  x  such  that  0  £  x^  _<  I 


(6.39) 


Mx)  *  V-  (x)  +  V.  (x). 

i  \2 


Proposition  6.3  Let  sk  be  a  permutation  of  {l,2,...,n}  such  that 
cksk(l)  -  cksk(2)  -  *•*  -  cksk(n). 

a.  If  i^  and  i2  have  permutations  such  that  si  ^  3  s-  j  *  l,...,n 
then  (6.38)  is  a  valid  aggregation. 

b.  If  i,  i^  and  i^  have  permutations  such  that  s^. ^  *  s^  3  s- 
j  3  l,...,n  then  (6.38)  is  a  valid  disaggregation. 

Proof  a.  If  s^  qj  3  s^  for  all  j  and  rows  i^  and  i^  are  aggre¬ 
gated  by  (6.38),  then  s^^  3  si  yj  *  j  3  1 . n.  Thus 

the  greedy  solution  to  ( 6 .24) -(6-27)  implies  that  (6.39)  holds, 
b.  The  proof  is  similar  to  that  of  a.  and  is  left  as  an  exercise.  Q 


The  following  proposition  shows  how  a  row  can  be  disaggregated. 


Proposition  6.4  Suppose  cis^_> 
Then  for  j  3  l,...,p  -  1 


■  -  ciSj(p-l)  >  c1si(p)  -  •”  -  cisj(n)' 


m 


•*,*4  U)  ‘  IS,  (p) 


*l2si2(j)  *  Clsi(j)  *  clst(p) 


and  for  j  •  p,...,n 


C11si  (j)  *  c1sf(j)  and  (j)  *  0 


is  a  valid  disaggregation  of  row  1. 

Proof  The  condition  of  Proposition  6.3b  holds.  Q 

We  can  apply  Proposition  6.4  recursively  to  disaggregate  a  row  i  into 

at  most  n  rows,  1j,...,1n  with  the  following  properties: 

(i)  c,  •  €  (0,r.  }  for  t  »  l,...,n,  (ii)  r.  >0  for  t  >  1,  (iii) 
tJ  Tt  h 

c,  .  *  r.  for  j  *  l,...,n,  (1v)  c,  ,  »  0  implies  c,  ,  *  0  for 

V  h  V  Vr3 

£  *  2,...,n-l. 

To  do  this,  suppose  c^  ^  *  ...  *  c^$  for  some  q  <  n. 

Apply  Proposition  6.4  with  p  *  q.  This  yields 


ev,  U>  '  cls,(<i)' 


j  *  1 , » . . »n 


Cts,0)  '  c(»,(q)  J,t . "  ‘  1 


J  *  q». « • *n« 


Row  i,  is  in  the  desired  form  and  if  c,„  »  c.„  j  *  l,...,q 

1  isi  (j )  iSjU+1) 

so  is  row  i».  Otherwise  let  i  be  the  largest  value  of  j  for  which 


^  ».  «/  *V  V  \" 


1MUCU\ 


m 


cis  (i-1)  >  c1s  (*)*  Now  app^y  Proposition  6.4  with  p  »  i  to  disaggregate 
row  ig.  Here  we  refer  to  the  two  new  rows  as  i2  and  i^.  Thus 


cVi2«> 


C1s1(Jl)  ‘  C1s,(q) 


j  *  1 . .  -  1 


J  *  q....,n 


V,3u> 


Cis1(j)  "  Cisi(t)  j  3  1 . l  -  1 


j  a  )  •  «  «  fHi 


Row  i2  is  now  in  the  desired  form  and  we  now  disaggregate  row  i3  if 
necessary.  Since  at  each  step,  the  row  to  be  disaggregated  has  at  least  one 
more  zero  than  the  previous  row,  the  procedure  takes  at  most  n  -  1  steps 
and  yields  at  most  n  rows  each  having  the  desired  property. 

Consider 


(ci l  ci2  C1 3  ci4^  *  (4  2  2  _1)* 


We  obtain 


<cijl  <^2  ‘^3  ■  '-1  -1  -1  -11 


<ci2l  c122  c123  ct24>  1  <5  3  3 
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Disaggregating  row  12  yields  (3330)  and  (200  0).  Thus  a  client 
represented  by  the  row  (422  -1)  can  be  replaced  by  3  equivalents  clients: 
(-1  -1  -1  -1),  (3330)  and  (200  0). 

Suppose  each  of  the  clients  1  e  I  is  disaggregated  in  this  way.  Now  we 


may  obtain  pairs  of  clients  say  1^  and  k^  such  that  j  *  0  and 

t 

c  .  *  0  if  and  only  if  j  «  T.  By  Proposition  6.3.  a.,  these  two  clients 
P*J 

can  be  aggregated  into  a  single  client  with  profit  c.  .  +  c  .  for  j  €  T 
and  0  profit  for  j  7  T.  Finally,  a  client  whose  profit  is  constant  for  all 


j  t  J  can  be  eliminated  from  the  problem  since  this  client  produces  the  same 
profit  for  all  feasible  x. 

To  summarize  this  discussion,  we  can  transform  the  matrix  C  Into  an 
equivalent  canonical  matrix  R  containing  at  most  min(m(n-l),2n-2)  rows. 
Each  row  of  R  represents  a  set  T  r  J\$  and  there  Is  a  profit  rT  for  all 

j  e  T  and  a  profit  of  zero  for  j  ^  T. 

In  our  example. 


W'H'-.'  O  K”  K"  V.W.'V.-.v'v  7^7  .1.  ■-.  -  „  . .  -_.  v  .  ..  ,  .-,  -  ,  -  . 
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We  will  now  use  this  transformation  to  obtain  a  reformulation  of  SIP. 
Given  x,  0  jC  x.  1  all  j  e  J,  the  problem  for  the  client  that  represents 

J 

the  set  T  that  is  equivalent  to  (6.24)-(6.27)  is 

VT(x)  -  rT(min(  7  x.,  I))  »  rT  -  rT(l  -  7  x.)\ 

'  1  j«T  J  r  r  j«T  J 

Let  ir-j.  be  the  fraction  of  client  T  not  served.  Then 

VT(x)  *  r-j.  +  max(-rTitT) 

itT  >  1  -  7  x. 

T”  jST  J 

■tij  £  0. 


Let 

be  the  collection  of  subsets  of 

J  that  are  represented  in  the  profit 

matrix 

R.  Then 

(6.40) 

Zlp  '  t Sr  fT 

+  max(-  r 

TU- 

le  / 

T"T  -  jj  Vj> 

(6.41) 

+  7  x  •  >  1 
JeT 

all  T  e  T 

(6.42) 

nT 

>0,  x.  >  0 

J 

all  T  £  T,  all  j  {  J. 

We  rewrite  the  objective  function  as 


*  T»  TW 


•«  * 


-I'  V 


(6.43) 


ZLP  *  ZLP  *  I  rr  3  "  min(  I  rT«T  +  I  fixi) 

LK  LK  T eT  1  TeT  T  1  j«J  J  J 


The  dual  of  (6.41)-(6.43)  is 


(6.44) 


Z*p  -  -  max  l  uT 
T IT  ' 


(6.45) 


?  uT  <  f  . 
Taj  T"  J 


all  j  €  J 


(6.46) 


0  <  Uj  <_  Pj  all  T  €  T 


The  linear  program  (6.44)-(6.46)  has  at  most  m(n-l)  variables  and  n 


constraints  plus  upper  bounds  on  the  variables.  It  is  the  most  compact 


linear  programming  formulation  we  know  and  has  the  structure  of  the  linear 


programming  relaxation  of  a  set  packing  problem.  In  our  example,  see  the 


matrix  R  given  above,  there  are  only  10  variables  and  6  constraints  other 


than  upper  bounds.  In  comparison,  the  original  linear  programming 


formulation  of  SLP  ( (6.1)— (6.4) )  has  32  variables  and  28  constraints.  In 


experimenting  with  some  k-median  problems,  Cornuejols,  Nemhauser  and  Wolsey 


(1980)  have  observed  that  when  the  simplex  method  (with  upper  bounds  treated 
implicitly)  is  applied  to  the  various  linear  programming  formulations,  the 


formulation  (6 . 44) -(6.4-6 )  was  by  far  the  best  one  in  terms  of  simplex  pivots 


and  running  time.  In  addition,  the  formulation  (6.44)-(6.46)  provides  a  nice 


interpretation  for  the  dual  descent  heuristic  given  in  Section  5.  We  leave 


V.V, 


this  as  an  exercise. 
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Summary 

This  section  has  emphasized  the  description  of  exact  algorithms  for  the 
solution  of  SIP.  These  algorithms,  and  effective  heuristics  for  the  dual, 
such  as  dual  descent,  can  be  viewed  as  subroutines  that  must  be  embedded 
within  a  branch-and-bound  code  based  on  enumeration  of  the  0-1  variables  x.. 

J 

A  second  viewpoint  is  based  on  a  specific  MIP  (mixed  integer  programming) 
formulation  that  can,  in  principle,  be  tackled  by  any  general  purpose  mixed 
integer  programming  package.  Here  candidates  are  the  original  formulation 
(1.1)— (1.4) ,  the  Benders  formulation  (6 . 34) -( 6. 36 )  with  x.  e  (0,1)  all 

J 

j  e  J,  and  the  new  canonical  formulation  (6.4Q)-(6.42)  with  x.  e  (0,1)  all 

J 

j  c  J.  Such  formulations  are  also  amenable  to  treatment  by  algorithms  other 
than  branch-and-bound— see  for  instance  the  cutting  planes  for  (1.1)-(1.4) 
developed  in  the  next  section.  However,  attempting  to  fit  any  of  these 
formulations  into  a  general  purpose  package  leads  to  difficulties  of  problem 
size  since  each  of  these  formulations  involves  G(mn)  constraints  and/or 
variables. 

Having  been  inundated  with  various  heuristics  and  algorithms,  the  reader 
has  the  right  to  ask  for  a  recommendation  on  the  appropriate  method  to  use. 
Unfortunately  there  is  no  simple  answer. 

Among  the  special  purpose  branch-and-bound  algorithms  0UAL0C  appears  to 
be  the  best.  It  is  a  generally  available  easy  to  use  FORTRAN  program  that  is 
very  fast  on  most  problems. 
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If  one  Insists  on  solving  the  strong  linear  program  to  optimality,  which 

may  be  necessary  for  the  harder  problems,  solving  the  lagrangian  dual  by 

subgradient  optimization  provides  an  easy  to  program  and  relatively  fast 

approach.  If  the  simplex  algorithm  is  to  be  used,  the  linear  program 

(6.44)-(6.46)  has  a  significant  advantage  In  size,  and  limited  computational 

experience  suggests  that  it  is  best  in  terms  of  time  and  pivots. 

The  main  advantage  of  using  a  general  purpose  MIP  code  is  that 

complicating  constraints  create  no  difficulties,  whereas  a  special  purpose 

code  becomes  unusable.  In  addition,  work  is  in  progress  on  general  purpose 

codes  which  will  be  capable  of  working  with  the  compact  formulation  "  y  -  •  <. 

mx.  and  generating  violated  variable  upper  bound  constraints  y.  .  £x.  as 
J  f  J  J 

needed;  see  Martin  and  Schrage  (1982)  and  Van  Roy  and  Wolsey  (1983).  This 
should  permit  such  codes  to  handle  medium  sized  UL's  as  general  MIP's. 


.v.v 


V-  1 


Poivhedrai  Results 


Extreme  Points  of  SIP 


We  have  observed  that  the  strong  linear  programming  relaxation 


m  n 


max  I  l  ciiyij  *  l  Vi 

i =1  j«l  j=l  J  J 


subject  to 


(7.1) 


j-i  J 

0  <_  yi j  <_  Xj  £  1  i  «  l,...,m  and  j  *  l,...,n. 


very  often  has  integral  optimal  solutions.  Although  this  phenomenon  is  not 

well  understood,  some  properties  of  the  polytope  (7.1)  are  known.  When 

m  _<  2  or  n  _<  2,  it  has  been  shown  by  Muckendi  (1975),  Krarup  and  Pruzan 

(1983)  and  Cho  et  al.  (1983)  that  all  the  extreme  points  of  (7.1)  are 

integral.  In  fact,  the  constraint  matrix  is  totally  unimodular  in  that  case. 

However,  even  when  m  *  n  »  3  the  strong  linear  programming  relaxation  may 

/°  1  l\ 

have  fractional  optimal  solutions.  For  example,  when  C  =  1  0  1  land  f  .  *  1 

Vi  i  o/  J 

for  j  *  1,2,3,  we  have  remarked  previously  that  x.  *  1/2  for  j  *  1,2,3 
and  yi j  *  1/2  for  i  t  j,  0  for  i  *  j,  where  i,  j  =  1,2,3  is  the  umaue 

optimal  solution. 


V»'. *•  . *»  *•  »*»  .  *  .  •  **  .*•  .  .*■  a  *  -***  .N  A  .  «“•  , *' '-V 


The  fractional  extreme  points  of  (7.1)  are  completely  characterized  by 
the  next  theorem.  For  a  given  non-integral  solution  (x,y)  of  (7.1)  lot 

J,  1  (j  s  J:  0  <  x.  <  1}  and 

*  w 

3  (i  €  I:  y^j  a  0  or  x^  for  all  j  and  y^  is  fractional  for  some  j} . 

Let  a.j  *  1  if  y.j  >  Q  and  0  otherwise,  and  denote  by  A  the  J I ^ j  x  jj^ 
matrix  whose  elements  are  a^  for  i  £  1^,  j  £  J^. 

Theorem  7.1  [(Cornuejol s,  Fisher  and  Nemhauser  (1977a)].  A  fractional  solution 
(x,y)  of  (7.1)  is  an  extreme  point  of  (7.1)  if  and  only  if 

(i)  x.  =  max  y .  .  for  all  j  £  J, 

J  ul  iJ  i 

(ii)  for  each  i  c  I,  there  is  at  most  one  j  s  J  with  0  <  y . .  <  x, 

I  J  J 

(iii)  the  rank  of  A  equals  |J^|. 

The  3  conditions  of  this  theorem  are  easily  verified  for  the  example  given 
above. 

Since  the  polyhedron  defined  by  (7.1)  has  many  fractional  extreme  points, 
the  type  of  objective  function  that  is  optimized  over  this  polyhedron  must 
play  an  important  role  in  the  attainment  of  integral  optimal  solutions. 
Frequently,  C  *  (c^.)  is  defined  over  a  network  with  the  property  that 
c^j  =  -d . (t i j  +  dj )  decreases  as  a  function  of  i,  the  further  i  is  from 
j  in  the  network.  (I.e.  if  i'  is  on  the  shortest  path  from  i  to  j 
in  the  network  then  c^  .  <  c.,j).  The  influence  of  this  property  on  the 
solution  of  SLP  for  tree  networks  'will  be  discussed  in  Section  8.  For  more 
general  graphs,  it  is  an  interesting  open  question. 


When  C  *  { c ^ j )  is  a  general  0,1  matrix  and  fj  *  1  all  j  e  J,  UL 

is  the  problem  of  finding  the  minimtri  nunber  of  columns  that  cover  all  the 

rows  of  the  maxtrix  C.  (A  set  S  of  columns  covers  row  i  if  c- .  =  1 

^  J 

for  at  least  one  j  «  S.)  This  problem  is  known  as  the  set  covering  problem 
and  often  has  fractional  optimal  solutions. 

Some  valid  inequalities  for  UL  that  remove  fractional  extreme  points  of 
the  SIP  polytope  are  given  in  the  next  theorem. 

Theorem  7.2  Cho  et  al.  (1983).  Let  8  be  a  k  x  k  nonsingular  0,1  matri 
such  that  8’*  e  >_  0,  where  e  is  a  column  vector  of  ones.  Index  the  rows 
and  columns  of  8  by  1^  c_  I  and  c_  J  where  |l^j  =  jJ^|  =  k.  Then 


uU 


A 

JeJ. 


b • .y. • 


l 

j«jk 


ili- 


eT8-1ej 


is  a  valid  inequality  for  UL.  It  cuts  off  at  least  one  fractional  extreme 

point  of  the  polytope  (7.1)  if  eT  8"^  e  is  not  integral. 

/0  i  i  \ 

For  example,  if  8=101  I  we  generate  the  constraint  y  +  y  + 

\1  1  0  /  12  13 

y?l  +  y23  +  y2^  +  y22  “  X1  “  *2  "  x3  —  wl11ch  cuts  the  fractional 


extreme  point  given  above. 

Conversely,  it  is  easy  to  show  that  the  family  of  valid  ineaualities 
defined  in  Theorem  7.2  cuts  off  all  the  fractional  extreme  points  of  the 
polytooe  (7.1).  However,  in  general  new  fractional  extreme  Doints  arise. 
We  now  turn  to  the  identification  of  valid  inequalities  for  UL  that 


define  facets  of  the  integer  polytooe. 


Facets  of  the  Integer  Polytope 

Let  Pm  n  be  the  polytope  defined  as  the  convex  hull  of  the  integer 
solutions  to  the  system  (7.1) 


r  n 


1  y. 

*  i 

i  *  l,...,m 

j»i  '' 

o  < 

'  VI 

< 1 

i  3  1 , . . .  ,m 

and 

j  s  1 

xr  yu 

<  (0,1) 

i  *  1 , . . .  ,m 

and 

j  -  1 

0  1  k 

A  set  of  k  +  1  points  w  ,  w  ,...,w  are  affinely  independent  if  the 

k  vectors  w*  -  w°,...,wk  -  w°  are  linearly  independent.  A  polytope  has 

dimension  k  if  it  contains  k  +  1  affinely  independent  points  but  not 

more.  An  affine  space  is  the  intersection  of  hyperplanes.  The  smallest 

affine  space  which  contains  a  polytope  P  is  called  its  affine  hul 1 .  The 

polytope  P  has  dimension  mn  +  n  -  m  and  affine  hull 
m,n 

n  mn 

{ (x,y)  e  R  x  R  :  >  y-  .  *  1  for  i  *  l,...,m). 

j3l 

A  face  of  a  polytope  P  is  a  set  F  =  P  n  {x:  ax  *  b)  where  ax  £  b 
is  satisfied  by  every  x  e  P  and  ax  <  b  for  at  least  one  x  e  P.  The 
inequality  ax  £  b  is  said  to  define  the  face  F.  Any  face  of  a  polytope  is 
itself  a  polytope.  When  its  dimension  is  one  less  than  that  of  the  polytope 
P,  the  face  F  is  called  a  facet.  To  describe  the  polytope  P  by  a  linear 
system,  it  suffices  to  have  a  description  of  its  affine  hull  and  one  defining 
inequality  for  each  facet  of  P. 


If  such  a  description  of  Pffl  n  were  know  then,  in  principle,  UL  could 

be  solved  as  a  linear  program  since  the  extreme  points  of  P  are 

m,n 

precisely  the  feasible  solutions  of  UL.  However,  a  complete  linear  system 
defining  Pffl  n  is  not  known  explicitly.  Even  if  one  were,  only  relevant 
portions  of  it  would  be  generated  to  solve  an  instance  of  UL,  i.e.,  the 
facets  of  P^  n  would  be  used  as  cutting  planes  in  the  spirit  of  Padberg  and 
Hong's  (1980)  and  Grotschel's  (1980)  work  on  the  traveling  salesman  problem. 
Whatever  the  algorithmic  use  of  a  partial  linear  description  of  P 

m,n 

the  first  step  is  to  Identify  some  of  its  facets. 

Theorem  7.3  The  following  inequalities  define  distinct  facets  of  P 
-  m,n 


(i) 

for  all 

1  e 

I.  j  e  J 

(ii) 

yij!° 

for  all 

i  c 

I.  j  €  J 

(iii) 

"J  A* 

for  all 

j  € 

J. 

These  facets  are  called  the  elementary  facets  of  Pm  n .  The  next  theorem 
provides  a  necessary  and  sufficient  condition  for  an  inequality  with  coeffi¬ 
cients  of  0  or  1  to  define  a  facet  of  P„  . 

m,n 

Assume  that  I'  £  I  and  J'  e_  J  are  two  nonempty  sets  and  that 
B  a  (bij)  '  £  I'«  j  e  J'  is  a  0,1  matrix  with  no  zero  row.  Consider  the 


inequality 


Id'  jeJ* 


ij-  ij 


jeJ*  J 


Oefine  the  graph  S  as  follows.  It  has  a  node  associated  with  each 
variable  y^j,  i  «  I,  j  t  J,  and  x^,  j  <  J.  We  will  use  the  same  notation 
for  a  node  and  Its  associated  variable.  For  all  i  e  I  and  j  e  J.  the 
node  y. .  is  joined  by  an  edge  to  the  node  x.  and  to  every  node  y..  for 

1  J  J  1 

k  *  j. 

Let  N'  be  the  set  of  nodes  {y^.J  1  «  1 1 ,  j  c  J*  :  {x^}  j  €  J'  and  let 
G*  be  the  subgraph  of  G  Induced  by  the  node  set  N'.  Given  a  graph  H  we 
denote  by  a(H)  the  maximum  size  of  a  stable  set  in  H  (a  stable  set  is  a 
set  of  mutually  nonadjacent  nodes).  Finally,  an  edge  e  of  H  Is  critical 
if  <z(H  -  e)  >  a(H) ,  where  H  -  e  denotes  the  graph  obtained  from  H  by 
removing  the  edge  e. 


Theorem  7.4  [Cornuejols  and  Thizy  (1982a)]  The  inequality  (7.2)  is  a  facet 


of  „  if  and  only  if  the  following  set  of  conditions  is  satisfied 
m,n  3 


(1)  r  -  a(G' )  -  |J| 

(ii)  G'  is  connected, 

(iii)  for  every  1  «  I',  j  «  J'  such  that  b^.  *  l,  the  edge  (x^.,  y^)  is 
critical , 

(iv)  for  every  j,  k  e  J',  there  exists  a  sequence  of  critical  edges 

(y*  4*  y,  ,  ).  (y<  ,  .  y<  .  )»*••. (y<  .  .  y,  ,  ).  (y,  ,  *  y<  k 

1J  ’n  1 2*1  '2*2  s-rs-2  ’s-rs-l  Vs-1  V 


(v)  for  every  i  e  I,  j  s  J  such  that  y  «  N',  the  inequality  a(G')  < 

ij 


a(G")  is  strict,  where  GH  denotes  the  subgraph  of  G  induced  by 


N'  j  {yi j}. 


These  necessary  and  sufficient  conditions  can  be  used  to  prove  the  next 


theorem,  which  provides  constructively  a  large  class  of  facets  for  the 
uncapacitated  plant  location  polytope. 

For  2£t  <  i  £  n,  define  B*j  ■  (b*j)  as  a  matrix  with  (*)  rows  and 
l  columns  whose  rows  consist  of  all  distinct  0,1  vectors  with  t  ones 
and  i  -  t  zeros. 

Theorem  7.5  [Cornuejols  and  Thizy  (1982a)]  For  any  pair  of  integers  i  and  t 
such  that  2  £  t  <  l  £  n  and  (*)  £  m,  and  any  sets  I*  £  I,  J'  <=_  j  such  that 
| I' |  a  (*)  and  Jj'|  =  i,  the  inequality 


n 


uv 


:  x.  <  (*) 

jcJ*  J  * 


i 


defines  a  facet  of  P_  . 

m,n 

23  /°  1  A 

For  example,  take  t  =»  2,  A  a  3,  B  *  1  0  1  )  ,  and 

\1  1  0  / 

{1,2,3}.  According  to  Theorem  7.5  we  get  the  facet 


I'  -  J*  - 


y12  +  y13  +  y21  +  y23  +  y31  +  y32  *  x!  “  x2  “  x3  -  1* 


which  is  Identical  to  the  valid  inequality  of  Theorem  7.2  that  we  obtained 
23 

above  with  B  »  B  . 

Other  facets  obtained  by  lifting  odd  holes  or  circulant  matrices  can  be 
found  in  [Cornuejols  and  Thizy  (1982a)].  Additional  material  can  be  found  in 
Cho,  Padberg  and  Rao  (1983). 


8y  a  special  case  of  UL,  we  mean  a  problem  of  the  form  (1.1)-(1.4)  all 
of  whose  instances  are  described  by  a  subfamily  of  objective  functions  (C,f). 
In  this  section,  we  consider  two  special  cases  of  UL  that  have  the  following 


significant  properties. 

1.  SLP  always  has  an  integral  optimal  solution 

2.  The  problem  can  be  solved  in  polynomial-time. 


A.  Economic  Lot  Sizing 

There  is  a  demand  di  in  period  i,  i  »  l,...,n.  The  fixed  cost  of 
producing  in  period  j  is  f.  ^  0.  The  variable  production  cost  is  p.. 
The  variable  storage  and  backorder  costs  are  ct  _>  0  and  cT  _>  0  respec- 

J  J 

tively.  In  UL,  y^  now  represents  the  fraction  of  the  demand  of  period  i 
produced  in  period  j,  and  x.  3  1  if  and  only  if  there  is  production  in 

J 

period  j. 


and 


eu  *  -<»j  *  ♦•••♦  cT-i>di  1f  1 
cij  •  -(pj  *  cI  *•••*  (f  1 


B.  The  Tree  Location  Problem 

Let  G  *  (V,E)  be  a  graph  with  node  set  V  and  edge  set  E  and  suDpose 
that  G  is  a  tree,  i.e.  there  is  a  unique  path  in  G  between  each  pair  of 
nodes  (G  is  connected  and  acyclic).  Here  the  nodes  represent  both  clients 
and  facilities.  The  cost  of  opening  the  jth  facility  is  f,  _>  0  all  j  e  V. 

J 

Associated  with  each  edge  e  e  E,  there  is  a  given  non-negative  distance.  The 


Kv: 


distance  between  nodes  1  and  j  is  d^.  *  sum  of  the  edge  distances  along 

the  unique  path  between  i  and  j,  for  all  i,  j  e  V.  There  is  also  a 

non-negative  scaling  function  w.  associated  with  each  node  1.  Let  c. .  * 

^  J 

j ,  i ,  j  e  V,  i  *  j  and  c^  »  0  all  i  e  V. 

The  induced  subgraph  of  G  generated  by  V.  V  is  the  graph  G.  = 
(V.,E.)  where  E.  *  {e  «  E:  both  end  nodes  of  e  are  in  V.}.  G.  is  said 

WWW  J  J 

to  be  a  subtree  of  G  If  G.  itself  is  a  tree. 

J 

Theorem  8.1  [Kolen  (1982)]  There  is  an  optimal  solution  to  the  tree  location 
problem  in  which  the  set  of  open  facilities  S  £  V  is  such  that  for  each 

j  «  S,  V.  *  {i  £  V:  node  1  is  served  by  node  j}  induces  a  subtree. 

J 

Moreover,  this  solution  is  also  optimal  to  the  linear  programming  relaxation 
of  the  tree  network  problem. 

A  similar  result  applies  to  the  lot  sizing  problem.  Consider  the  tree 
G  =  (V,E)  where  V  *  {1,2,... ,n}  and  E  *  {(i,i+l):  i  *  l,...,n-l>.  Here 
G  is  simply  a  path  from  node  1  to  node  n  so  that  V'  =  V  is  a  subtree 

or  a  path  if  and  only  if  V'  *  {1,  1+l,...,k}  for  some  1  and  k,  1  £i  £n 

and  k  _>  i . 

Theorem  8.2  [Krarup  and  Bilde  (1977)1  There  is  an  optimal  solution  to  the  lot 
sizing  problem  in  which  the  set  of  periods  having  positive  production  Sc  V 
is  such  that  for  each  j  s  S,  V.  »  {i  €  V:  period  i  is  served  by  production 

w 

in  period  j}  induces  a  path.  Moreover,  this  solution  is  also  optimal  to  the 
linear  programming  relaxation  of  the  tree  network  problem. 

The  fact  that  an  optimal  solution  to  these  problems  induces  subtrees  that 
partition  V  is  not  surprising.  In  the  tree  location  problem,  suppose  that 


ki: 


IP 


t »  .!*  »  TmM 


node  1  Is  on  the  path  joining  nodes  j  and  k  and  node  j  serves  node  k 
but  not  node  1.  Suppose  node  1  Is  served  by  node  j'  *  j.  Then  the  cost 
from  this  part  of  the  solution  is  d^k  +  fj  +  dj'i  +  fj'*  instead,  node  1 

Is  served  by  node  j,  the  cost  Is  d^  +  djk  +  fj  +  fj'* 

Thus  If  d^  £  d j , ^ ,  Me  can  serve  1  and  all  nodes  after  i  that  are 
being  served  by  j'  without  increasing  the  cost.  This  reduces  the  number 
of  subtree  violations  by  one.  Otherwise  d..  >  d.,.,  which  implies 

J  i  J  i 


d.,  =  d . .  +  d. ,  >  d . - .  +  d. ,  *  d . . ,  . 

jk  ji  lk  j'i  ik  j'k 


This  inequality  implies  that  the  solution  in  which  j'  serves  k  and  only 
nodes  after  k  that  are  currently  being  served  by  j  costs  less  than  the 
original.  It  also  reduces  the  number  of  subtree  violations  by  at  least  one. 

A  similar  argument  proves  the  result  for  the  economic  lot  sizing  problem. 
Both  of  these  results  suggest  that  the  ability  to  partition  the  solution 
into  subtrees  is  crucial  and  leads  us  to  consider  the  following  generalization. 

C.  The  Tree  Partitioning  Problem  Given  a  tree  graph  G  *  (V,E)  and  a  node  by 
node  matrix  with  elements  all  i,  j  €  V,  let  the  weight  of  a  subtree 

G.  *  (V.,E.)  be  w(G.)  *  max  (  "  y.J.  Find  a  partition  of  G  into  sub- 

J  J  J  J  |<€y  .  •)  “y  _  1  * 

J  J 

trees  such  that  the  sum  of  the  weights  over  all  subtrees  in  the  solution  is 


maximum. 


To  model  the  lot  sizing  problem  as  a  tree  partitioning  problem  we  take 


*  -(fj  +  pjdj> 

tu  •  *  cj  *•••*  ci.i>dt  ,f  1  >  j 

Yij  »  -CPj  *  Cj  C-+1)d,  If  i  <  j. 

To  model  the  tree  location  problem  as  a  tree  partitioning  problem  we  take 

TjJ.  -fj  and  Tlj  *  — j  if  1  *  J- 

We  now  formulate  the  tree  partitioning  problem  as  an  integer  program  in  a 

manner  that  establishes  Its  connection  to  UL.  If  j*  e  V.  is  such  that 

J 

arg  max  (  "  y..  )  *  j*.  we  say  that  j*  is  the  root  of  subtree  G.»  Let 

k«V:vj  J 

■  1  if  i  €  V  is  In  a  subtree  rooted  at  j  and  y.^  *  0  otherwise. 
Then  the  tree  partitioning  problem  can  be  formulated  as 

<8-D  "»*  I  I 

i  j 

(8.2)  :  *  1  all  1  e  V 

j  J 

(8.3)  y. ,  -  y...  <0  all  i,  i',  j  £  V  such  that  i' 

J  J  precedes  i  on  a  path  from  j  to  i 

(8.4)  y..  e  {0,1}  all  1,  j  £  V. 


Constraint  (8.3)  guarantees  that  If  j  is  the  root  of  a  tree  that  contains 
i  then  the  tree  must  also  contain  1*.  Constraint  (8.2)  guarantees  that  each 


node  must  be  in  exactly  one  tree. 

The  linear  programming  relaxation  of  this  Integer  program  Is  obtained  by 
replacing  (8.4)  by 


(8.5) 


y-i  j  > 0  a11 


Theorem  8.3  [Barany,  Edmonds  and  Wolsey  (1983)]  The  polyhedron  defined  by 
(8.2),  (8.3)  (8.5)  has  only  Integral  extreme  points.  Hence  for  any  objective 
function  (8.1)  the  solution  to  the  linear  programming  relaxation  is  integral. 

This  model  for  the  tree  partitioning  problem  resembles  the  model  (1.1)- 
(1.4)  for  UL  if  we  think  of  the  y^ 's  as  x^'s  and  replace  (8.3)  by  (1.3). 
In  fact,  we  can  represent  solutions  of  UL  as  a  collection  of  subtrees  that 
partition  a  graph  G,  but  unfortunately  G  is  not  a  tree.  In  the  graph  of 
Figure  8.1,  V1  represents  the  set  of  clients  and  V2  the  set  of  facilities. 


Figure  8.1 


A  solution  to  the  problem  is  a  set  of  subtrees,  each  of  which  is  rooted  at 
a  node  in  V2.  The  edges  from  the  root  of  a  tree  to  nodes  in  (solid 
edges)  show  the  clients  being  served  by  the  facility  that  corresponds  to  the 
root.  A  node  in  V2  that  is  not  a  root  corresponds  to  an  unused  facility 
(dashed  edges).  We  set  y.  ^  =  c^  for  i  €  and  j  *  V^.  To  eliminate 
the  possibility  of  nodes  in  being  roots  we  set  y^  *  -M  for  j  £  V1 

(M  is  a  large  positive  number)  and  to  represent  the  fixed  costs  we  set 
y  »  -f  for  j  £  V?.  Finally  to  accommodate  unused  facilities  we  set 
Yjk  *  0  if  j  and  k  are  in  V2. 

i  2i 

We  close  this  section  by  giving  an  0(  V  )  dynamic  programming 
algorithm  for  solving  the  tree  partitioning  problem.  Another  very  general 
dynamic  programming  algorithm  for  location  problems  on  trees  can  be  found  in 
Megiddo,  Zemel  and  Hakimi  (1983). 


Given  the  tree  G  s  (V,E)  we  choose  an  arbitrary  root  r.  This  induces 
a  partial  order  on  V.  For  all  v  e  V,  let  V(v)  =  {t:  v  is  on  the  unique 
path  between  r  and  t}  and  S(v)  »  { t :  v  is  the  node  that  precedes  t  on 
the  unique  path  between  r  and  t}.  Let  be  the  tree  induced  by  V(v) 

and  g(v)  =  optimal  weight  partition  for  the  tree  T  . 

The  idea  of  the  algorithm  is  to  calculate  g(r)  recursively  by  deter- 
ming  g(v)  from  the  g(w)  for  all  w  e  S(v),  A  node  w  is  said  to  be  a 

leaf  of  T  if  S(w)  =  0.  Note  that  we  can  begin  the  recursion  with  g(w)  = 

rww  ^or  leaves.  Before  giving  the  general  recursion  equation,  we  need 

one  more  definition.  Let  ( v )  *  optimal  weight  partition  of  T^  when  v 

is  served  by  node  u  and  u  may  not  necessarily  be  an  element  of  V(v). 

Then 

(8.6)  g ( v )  »  max  g(v) 

ueV(v)  u 

and  if  w  is  a  leaf  of  T 

r 

(8.7)  gu(w)  =  ywu  f°r  all  u  e  V. 

Now  suppose  we  are  given  9u(w)  all  w  e  S(v)  and  all  u  «•  V.  The  calcu¬ 
lation  of  g  (v)  divides  into  two  cases  as  shown  in  Figure  3.2 


Figure  8.2 


If  u  ^  V(v)  or  u  *  v  and  u  serves  v,  then  u  will  also  serve 
w  e  S( v )  if  qu(w)  _>  g(w).  Hence 

(3-8)  9,.(v)  =  y  +  l  max{g  (w)  ,g(w)> . 

U  VU  wcS(v)  u 

If  u  c  V(v)\{v},  then  u  e  V(w*)  for  some  w*  e  S(v).  Hence  if  u 
v  then  u  serves  w*.  Thus 

(3.9)  g  (v)  =  y  +  l  max{g  (w),g(w)}  +  glw*'. 

u  VU  weS(v)\{w*}  u  u 

An  Example 

We  consider  a  tree  location  problem  on  the  qraoh  of  Figure  3.3. 
numbers  on  the  edges  are  the 
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Hence 

we  obtain  the  matrix 

A7 

-4 

-3 

-7 

-e\ 

/  "8 

-4 

-14 

-6 

‘4  \ 

r  -  -6 

-14 

-11 

-20 

-18 

V  -7 

-3 

-10 

-6 

-5  j 

\-l2 

-4 

-18 

-10 

-8  / 

For  w  =  3,4,5  the  bottom  3  rows  of  the  matrix  y  give  9u(w),  see 
For  node  2,  (8.8)  and  (8.9)  yield 

g  1  ( 2)  *  y21  +  »nax{g1(4),g(4)}  +  max{g1(5),g(5)} 

-8  +  max(-7,-6)  +  max(-12,-8)  =  -22 
g^(2)  3  -4  +  max(-3,-6)  +  max(-4,-8)  =  -11 

g3 (2)  »  -14  +  max(-10,-6)  +  max(-18,-8)  3  -28 

g4(2)  3  y24  ♦  max{g4(5),g(5)}  +  g(4) 

*  -6  +  max{ -10,-3}  +  (-6)  =  -20 

g j. ( 2)  3  -4  +  max (-5, -6}  -  (-8)  3  -17. 


9(2)  »  max{g2(2),g4(2),g5(2)}  =  max{-ll,-20,-17}  *  g2(2)  =  -11. 


For  node  1  we  obtain 


9j(l)  * 
92(D  - 

g3d)  ■ 
g4(!)  3 

95(D  - 


-7  +  max{ -22,-11}  +  max{ -6,-11}  = 
-4  +  max{-14,-ll}  +  (-11)  =  -26 
-3  +  max{ -28,-11}  +  (-11)  =  -25 
-7  +  max{ -20,-11}  +  (-20)  *  -38 
-6  +  max{ -13,-11}  +  (-17)  *  -34  . 


Hence  g(l)  »  g^l)  =  -24,  where  gx(l)  »  yu  +  g(2)  +  g1(3) . 

Thus  node  1  serves  itself  and  node  3.  Since  g(2)  =  g2(2)  »  y22  * 
g2(4)  +  g2(5),  node  2  serves  itself  and  nodes  4  and  5. 


9.  Submodularity 

As  defined  at  the  outset  of  this  chapter,  UL  is  the  combinatorial 
problem 


max  z(S) 
SeJ 


where 


(9.1) 


z(S)  »  l  max  c..  -  m 
i»l  US  1J  jcS 


is  the  profit  made  when  the  set  S  of  facilities  is  open.  A  very  important 
property  of  the  set  function  z  is  its  submodularity,  A  function  w 
defined  on  the  subsets  of  a  finite  set  J  is  submodular  if 


w(Su{k})  -  w(S)  £  w(Ru{k})  -w(R)  for  all  k  ^  S  and  Rc  Sc  J  -  {k}. 


The  fact  that  the  profit  function  z  is  submodular  was  observed  by 
Spielberg  (1969a),  Babayev  (1974),  Frieze  (1974)  and  Fisher,  Nemhauser  and 
Wolsey  (1978).  It  means  that  the  additional  profit  that  can  be  made  by 
opening  a  facility  In  location  k  when  a  set  S  is  already  open  in  other 
locations  is  a  nonincreasing  function  of  S  with  respect  to  set  inclusion. 
The  larger  S,  the  smaller  the  profit  of  establishing  a  new  facility.  This 
is  proved  formally  in  the  next  theorem. 


Theorem  9.1  The  profit  function  z  is  submodular. 


Proof  Let  R  =  S  c  J  -  {k}.  For  all  i  ■ 


max  c,  .  -  max  c.,  *  max(0,  c..  -  max  c, .) 
j«Su{k}  j£S  1J  1IC  jcS  1J 


<  max(0,  c.,  *  max  c, A  ■  max  -  max  c . . 

K  jeR  J  jcRu{k}  jeR 


where  the  inequality  follows  from  max  c. .  >  max  c,  .. 

jeS 

By  summing  these  inequalities  for  all  i ,  we  obtain 


m  mm  m 

max  c..  -  l  max  cH  <  "  max  c .  .  -  "  max  c... 

1-1  j«Su{k}  1*1  j€S  “  i-1  jeRu{k}  i*l  j€R 


Hence 


z(Su{k} )  -  z(S)  _<  z(Ru{k})  -  z(R).  Q 


Thus  UL  is  a  special  case  of  the  more  general  problem 


(9.2)  max  (z(S):  z  submodular}. 

ScJ 

We  can  apply  the  greedy  and  interchange  heuristics  to  (9.2),  we  can  formulate 
(9.2)  as  an  integer  program  and  many  of  the  results  that  we  have  given  for  UL 
extend  to  (9.2)  and,  in  particular,  to  the  capacitated  location  problem  which 
is  another  special  case  of  it.  We  will  not  elaborate  on  these  results  here, 
but  refer  the  interested  reader  to  Fisher,  Nemhauser  and  Wolsey  (1978), 
Nemhauser  and  Wolsey  (1978),  Nemhauser,  Wolsey  and  Fisher  (1978),  Cornuejols, 
Nemhauser  and  Wolsey  (1980)  and  Nemhauser  and  Wolsey  (1981). 
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